Wavefront propagation from one plane to another

ABSTRACT

The present invention provides methods, systems and software for scaling optical aberration measurements of optical systems. In one embodiment, the present invention provides a method of reconstructing optical tissues of an eye. The method comprises transmitting an image through the optical tissues of the eye. Aberration data from the transmitted image is measured across the optical tissues of the eye at a first plane. A conversion algorithm is applied to the data, converting it to corrective optical power data that can be used as a basis for constructing a treatment for the eye at a second plane.

CROSS-REFERENCE TO RELATED APPLICATIONS

This application is a divisional of U.S. patent application Ser. No.13/023,871, filed Feb. 9, 2011, which is a continuation of U.S. patentapplication Ser. No. 12/421,246, filed Apr. 9, 2009, which is adivisional of U.S. patent application Ser. No. 11/736,353, filed Apr.17, 2007, which also claims the benefit of U.S. Patent Application No.60/826,636 filed Sep. 22, 2006, and which is a continuation-in-part ofU.S. patent application Ser. No. 11/032,469 filed Jan. 7, 2005, whichclaims the benefit of U.S. Patent No. 60/550,514 filed Mar. 3, 2004. Thefull disclosure of each of these filings is incorporated herein byreference.

BACKGROUND OF THE INVENTION

The present invention generally relates to scaling optical aberrationmeasurements of optical systems. More particularly, the inventionrelates to improved methods and systems for processing optical powermeasurements taken at a first plane and converting those powermeasurements to corrective optical power measurements that can be usedat a second plane. The present invention may be useful in any of avariety of ocular treatment modalities, including ablative laser eyesurgery, contact lenses, spectacles, intraocular lenses, and the like.

Known laser eye surgery procedures generally employ an ultraviolet orinfrared laser to remove a microscopic layer of stromal tissue from thecornea of the eye. The laser typically removes a selected shape of thecorneal tissue, often to correct refractive errors of the eye.Ultraviolet laser ablation results in photodecomposition of the cornealtissue, but generally does not cause significant thermal damage toadjacent and underlying tissues of the eye. The irradiated molecules arebroken into smaller volatile fragments photochemically, directlybreaking the intermolecular bonds.

Laser ablation procedures can remove the targeted stroma of the corneato change the cornea's contour for varying purposes, such as forcorrecting myopia, hyperopia, astigmatism, and the like. Control overthe distribution of ablation energy across the cornea may be provided bya variety of systems and methods, including the use of ablatable masks,fixed and moveable apertures, controlled scanning systems, eye movementtracking mechanisms, and the like. In known systems, the laser beamoften comprises a series of discrete pulses of laser light energy, withthe total shape and amount of tissue removed being determined by theshape, size, location, and/or number of laser energy pulses impinging onthe cornea. A variety of algorithms may be used to calculate the patternof laser pulses used to reshape the cornea so as to correct a refractiveerror of the eye. Known systems make use of a variety of forms of lasersand/or laser energy to effect the correction, including infrared lasers,ultraviolet lasers, femtosecond lasers, wavelength multipliedsolid-state lasers, and the like. Alternative vision correctiontechniques make use of radial incisions in the cornea, intraocularlenses, removable corneal support structures, and the like.

Known corneal correction treatment methods have generally beensuccessful in correcting standard vision errors, such as myopia,hyperopia, astigmatism, and the like. However, as with all successes,still further improvements would be desirable. Toward that end,wavefront measurement systems are now available to accurately measurethe refractive characteristics of a particular patient's eye. Oneexemplary wavefront technology system is the VISX WaveScan® System,which uses a Hartmann-Shack wavefront lenslet array that can quantifyaberrations throughout the entire optical system of the patient's eye,including first- and second-order sphero-cylindrical errors, coma, andthird and fourth-order aberrations related to coma, astigmatism, andspherical aberrations.

Wavefront-driven vision correction has become a top choice for higherquality vision, after a series of significant development in theresearch of the wavefront technology (Liang, J. et al., J. Opt. Soc. Am.A 11:1949-1957 (1994); Liang, J. et al., J. Opt. Soc. Am. A 14:2873-2883(1997); Liang, J. et al., J. Opt. Soc. Am. A 14:2884-2892 (1997);Roorda, A. et al., Nature 397:520-522 (1999)). Although the ocularaberrations can be accurately captured, several factors need to beconsidered when they are corrected using, say, the refractive surgicaltechnique. The first of such factors is the relative geometrictransformation between the ocular map when the eye is examined and theocular map when the eye is ready for laser ablation. Not only can theeye have x- and y-shift between the two maps, but it can also havepossible cyclo-rotations (Walsh, G. Ophthal. Physiol. Opt. 8:178-182(1988); Wilson, M. A. et al., Optom. Vis. Sci. 69:129-136 (1992);Donnenfeld, E. J. Refract. Surg. 20:593-596 (2004) Chernyak, D. A. J.Cataract. Refract. Surg. 30:633-638 (2004)). Such problems have beenstudied by Guirao et al. (Guirao, A. et al., J. Opt. Soc. Am. A18:1003-1015 (2001)). Another problem deals with the pupil size change(Goldberg, K. A. et al., J. Opt. Soc. Am. A 18:2146-2152 (2001);Schwiegerling, J. J. Opt. Soc. Am. A 19:1937-1945 (2002); Campbell, C.E. J. Opt. Soc. Am. A 20:209-217 (2003)) using Zernike representation(Noll, R. J. J. Opt. Soc. Am. 66:203-211 (1976); Born, M. et al.,Principles of Optics, 7th ed. (Cambridge University Press, 1999)).Because of the analytical nature and the popularity of Zernikepolynomials, this problem has inspired an active research recently (Dai,G.-m. J. Opt. Soc. Am. A 23:539-543 (2006); Shu, H. et al., J. Opt. Soc.Am. A 23:1960-1968 (2006); Janssen, A. J. E. M. et al., J. Microlith.,Microfab., Microsyst. 5:030501 (2006); Bará, S. et al., J. Opt. Soc. Am.A 23:2061-2066 (2006); Lundström, L. et al., J. Opt. Soc. Am. A(accepted)).

Wavefront measurement of the eye may be used to create a high orderaberration map or wavefront elevation map that permits assessment ofaberrations throughout the optical pathway of the eye, e.g., bothinternal aberrations and aberrations on the corneal surface. Theaberration map may then be used to compute a custom ablation pattern forallowing a surgical laser system to correct the complex aberrations inand on the patient's eye. Known methods for calculation of a customizedablation pattern using wavefront sensor data generally involvesmathematically modeling an optical surface of the eye using expansionseries techniques. More specifically, Zernike polynomials have beenemployed to model the optical surface, as proposed by Liang et al., inObjective Measurement of Wave Aberrations of the Human Eye with the Useof a Hartmann-Shack Wave-front Sensor, Journal Optical Society ofAmerica, July 1994, vol. 11, No. 7, pp. 1949-1957, the entire contentsof which is hereby incorporated by reference. Coefficients of theZernike polynomials are derived through known fitting techniques, andthe refractive correction procedure is then determined using the shapeof the optical surface of the eye, as indicated by the mathematicalseries expansion model.

There is yet another problem that remains unaddressed. Opticalmeasurements such as wavefront measurements are often taken at ameasurement plane, whereas optical treatments may be needed at atreatment plane that is different from the measurement plane. Thus,power adjustments are often used when devising optical treatments forpatients. For example, power adjustments can be used by optometristswhen prescribing spectacles for patients. Typically, refractivemeasurements are made by an optometer at a measurement plane somedistance anterior to the eye, and this distance may not coincide withthe spectacle plane. Thus, the measured power corresponding to themeasurement plane may need to be converted to a corrective powercorresponding to the spectacle or treatment plane. Similarly, whenwavefront measurements are obtained with wavefront devices, in manycases the measured map is conjugated to the pupil plane, which is notthe same as the corneal plane or spectacle plane. To enhance theeffectiveness of a refractive surgical procedure, vertex correction maybe needed to adjust the power of the measured maps. Yet there remains alack of efficient methods and systems for such power conversions. Inother words, when the ocular aberrations are captured, they are often onthe exit pupil plane. However, when the correction is applied, it isoften on a different plane. For example, for refractive surgery, it ison the corneal plane. For contact lens, it is on the anterior surface ofthe contact lens. For intraocular lens, it is on the lens plane. And forspectacles, it is on the spectacle plane. Traditionally, for low orderspherocylindrical error, a vertex correction formula can be applied(Harris, W. F. Optom. Vis. Sci. 73:606-612 (1996); Thibos, L. N. S. Afr.Optom. 62:111-113 (2003)), for example, to archive the power correctionfor the so-called conventional treatment for refractive surgery. Thesame formula can be applied to the power calculation for visioncorrection using the contact lens, intraocular lens, and spectacles.However such formulas may not be useful in some cases, for example wherethere are high order ocular aberrations to be corrected. Hence, newformulas are needed to represent the ocular aberrations when they arepropagated to a new plane.

Therefore, in light of above, it would be desirable to provide improvedmethods and systems for processing optical data taken at a measurementplane and converting that optical data to corrective optical data thatcan be used at a treatment plane.

BRIEF SUMMARY OF THE INVENTION

The present invention provides methods and systems for processingoptical power measurements taken at a first plane and converting thosepower measurements to corrective optical power measurements that can beused at a second plane.

In wavefront-driven vision correction, although ocular aberrations areoften measured on the exit pupil plane, the correction is applied on adifferent surface. Thus there is a need for new systems and methods thataccount for the changes occurring in the propagated wavefront betweentwo surfaces in vision correction. Advantageously, embodiments of thepresent invention provide techniques, based on geometrical optics andZernike polynomials for example, to characterize wavefront propagationfrom one plane to another. In some embodiments, properties such as theboundary and the magnitude of the wavefront can change after thepropagation. The propagation of the wavefront can be nonlinear. Taylormonomials can be effectively used to realize the propagation. Approachesused to identify propagation of low order aberrations can be verifiedwith a classical vertex correction formula. Approaches used to identifypropagation of high order aberrations can be verified with Zemax®. Thesetechniques can be used with the propagation of common opticalaberrations, for example. Advantageously, the techniques describedherein provide improved solutions for wavefront driven vision correctionby refractive surgery, contact lens, intraocular lens, and spectacles.

Embodiments encompass techniques for treating an ocular wavefront whenit is propagated from one plane or surface to another. Zernikepolynomials can be used to represent the ocular wavefront; they areorthonormal over circular pupils (Born, M. et al., Principles of Optics,7th ed. (Cambridge University Press, 1999)). Taylor monomials can beused for the calculation of the wavefront slopes (Riera, P. R. et al.,Proc. SPIE 4769, R. C. Juergens, ed., 130-144 (2002); Dai, G.-m. J. Opt.Soc. Am. A 23:1657-1666 (2006); Dai, G.-m. J. Opt. Soc. Am. A23:2970-2971 (2006)). In some embodiments, Taylor monomials can be usedfor wavefront propagation. Before and after the propagation, Zernikepolynomials can be converted to and from Taylor monomials usingavailable conversion formulas (Dai, G.-m. J. Opt. Soc. Am. A23:1657-1666 (2006)). Some embodiments encompass the use of an orderingconvention for Zernike polynomials such as the ANSI standard (AmericanNational Standard Institute, Methods for reporting optical aberrationsof eyes, ANSI 280.28-2004 (Optical Laboratories Association, 2004),Annex B, pp. 1928). It has been discovered that in some embodiments,high order aberrations may undergo certain changes as a result ofwavefront propagation. For example, high order aberrations such as comamay present an elliptical, a bi-elliptical, a four-fold elliptical, oranother noncircular shape or boundary. Embodiments provide solutions tothe current needs for wavefront propagation techniques. Embodimentsinclude approaches for addressing low and high order ocular aberrationsas they propagate. Embodiments also encompass verification techniquesfor low order and high order aberration approaches. Further, embodimentsinclude the propagation of wavefronts with single-term aberrations.

In a first aspect, embodiments of the present invention provide a methodof calculating a refractive treatment shape for ameliorating a visioncondition in an eye of a patient. The method can include, for example,determining a measurement surface aberration corresponding to ameasurement surface of the eye, where the measurement surface aberrationincludes a measurement surface boundary and a measurement surfacemagnitude. The method can also involve determining a propagationdistance between the measurement surface of the eye and a treatmentsurface, and determining a treatment surface aberration based on themeasurement surface aberration and the propagation distance. Thetreatment surface aberration can include a treatment surface boundaryand a treatment surface magnitude. The method can also includecalculating the refractive treatment shape based on the treatmentsurface aberration. In some cases, the refractive treatment shape isconfigured to ameliorate a high order aberration of the measurementsurface aberration. In some cases, a difference between the treatmentsurface magnitude and the measurement surface magnitude is proportionalto the propagation distance. In some cases, a difference between thetreatment surface magnitude and the measurement surface magnitude isproportional to a direction factor. In some cases, a difference betweenthe treatment surface magnitude and the measurement surface magnitude isinversely proportional to a dimension of the measurement surfaceboundary. In some cases, a difference between the treatment surfacemagnitude and the measurement surface magnitude is inverselyproportional to a squared radius of the measurement surface boundary.The measurement surface of the eye can correspond to a pupil plane ofthe eye, and the treatment surface can correspond to a corneal plane ora spectacle plane of the eye. In some aspects, the measurement surfaceaberration includes a wavefront measurement surface aberration, and thetreatment surface aberration includes a wavefront treatment surfaceaberration. In some aspects, the measurement surface boundary includes awavefront measurement surface boundary, and the treatment surfaceboundary includes a wavefront treatment surface boundary. A measurementsurface magnitude may include a set of measurement surface coefficients,and a treatment surface magnitude may include a set of treatment surfacecoefficients. Optionally, a measurement surface magnitude may include aset of measurement surface wavefront coefficients, and a treatmentsurface magnitude may include a set of treatment surface wavefrontcoefficients.

In another aspect, embodiments of the present invention encompassmethods of calculating a refractive treatment shape for ameliorating avision condition in an eye of a patient. A method may include, forexample, determining a first wavefront measurement corresponding to thepupil plane of the eye, where the first wavefront measurement includes afirst wavefront boundary and a first set of wavefront coefficients. Themethod may also include determining a propagation distance between thepupil plane of the eye and a treatment surface, and determining apropagated wavefront measurement corresponding to the treatment surfacebased on the first wavefront measurement and the propagation distance.The propagated wavefront measurement can include a second wavefrontboundary and a second set of wavefront coefficients. The method may alsoinclude calculating the refractive treatment shape based on thepropagated wavefront measurement. In some cases, the treatment surfacecorresponds to a corneal surface, a spectacle surface, a scleral lenssurface, a contact lens surface, or an intraocular lens surface. Methodsmay also involve applying the refractive treatment shape to the eye ofthe patient to ameliorate the vision condition. In some cases, the erefractive treatment shape is applied to the eye of the patient in aselected treatment modality. The example, the method can encompassablating a corneal surface of the eye to provide a corneal surface shapethat corresponds to the refractive treatment shape, providing thepatient with a contact lens that has a shape that corresponds to therefractive treatment shape, providing the patient with a spectacle thathas a shape that corresponds to the refractive treatment shape,providing the patient with a scleral lens that has a shape thatcorresponds to the refractive treatment shape, or providing the patientwith an intraocular lens that has a shape that corresponds to therefractive treatment shape.

In some aspects, embodiments of the present invention include systemsfor generating a refractive treatment shape for ameliorating a visioncondition in an eye of a patient. A system may include, for example, aninput module that accepts a measurement surface aberration correspondingto a measurement surface of the eye, where the measurement surfaceaberration includes a measurement surface boundary and a measurementsurface magnitude. The system can also include a transformation modulethat derives a treatment surface aberration corresponding to a treatmentsurface of the eye. The treatment surface aberration may be based on themeasurement surface aberration and a propagation distance between themeasurement surface and a treatment surface. The treatment surfaceaberration may include a treatment surface boundary and a treatmentsurface magnitude. A system may also include an output module thatgenerates the refractive treatment shape based on the treatment surfaceaberration. In some system embodiments, a difference between thetreatment surface magnitude and the measurement surface magnitude isproportional to the propagation distance. In some embodiments, adifference between the treatment surface magnitude and the measurementsurface magnitude is proportional to a direction factor. In someembodiments, a difference between the treatment surface magnitude andthe measurement surface magnitude is inversely proportional to adimension of the measurement surface boundary. In some embodiments, adifference between the treatment surface magnitude and the measurementsurface magnitude is inversely proportional to a squared radius of themeasurement surface boundary. A treatment surface can correspond to acorneal surface, a spectacle surface, a scleral lens surface, a contactlens surface, or an intraocular lens surface, for example.

In a further aspect, embodiments of the present invention encompassmethods for characterizing an electromagnetic field that is propagatedfrom a first surface to a second surface, and systems for carrying outsuch methods. Exemplary methods may involve determining a first surfacecharacterization of the electromagnetic field corresponding to the firstsurface, where the first surface characterization includes a firstsurface field strength. Methods may also involve determining apropagation distance between the first surface and a second surface, anddetermining a second surface characterization of the electromagneticfield based on the first surface characterization and the propagationdistance. In some cases, the second surface characterization includes asecond surface field strength. In some cases, the first surface fieldstrength includes a first surface field phase, and the second surfacefield strength includes a second surface field phase.

In one aspect, the present invention provides a method of determining arefractive treatment shape for ameliorating a vision condition in apatient. The method comprises measuring a wavefront aberration of an eyeof the patient in order to provide a measurement surface aberration,deriving a treatment surface aberration of the eye based on themeasurement surface aberration, and determining the refractive treatmentshape based on the treatment surface aberration of the eye. Thewavefront aberration can correspond to a measurement surface that isdisposed at or near a pupil plane of the eye, and the treatment surfaceaberration can correspond to a treatment surface that is disposed at ornear an anterior surface of a cornea of the eye. The treatment surfaceaberration may be derived using a difference between the measurementsurface and the treatment surface.

In another aspect, the present invention provides a method ofameliorating a vision condition in a patient. The method comprisesmeasuring a wavefront aberration of an eye of the patient in order toprovide a measurement surface aberration, deriving a treatment surfaceaberration of the eye from the measurement surface aberration,determining a refractive treatment shape based on the treatment surfaceaberration of the eye, and applying the refractive treatment shape tothe eye of the patient to ameliorate the vision condition. The wavefrontaberration can correspond to a measurement surface that is disposed ator near a pupil plane of the eye. The treatment surface aberration cancorrespond to a treatment surface that is disposed at or near ananterior corneal surface of the eye, or a treatment surface thatcorresponds to a spectacle plane of the eye. Relatedly, the treatmentsurface may be disposed posterior to a pupil plane of the eye. Thetreatment surface aberration may be based on a difference between themeasurement surface and the treatment surface.

In a related aspect, the refractive treatment shape can be applied tothe eye of the patient in a variety of treatment modalities. Forexample, the treatment shape can be applied by ablating a cornealsurface of the patient to provide a corneal surface shape thatcorresponds to the refractive treatment shape. The treatment shape mayalso be applied by providing the patient with a contact lens that has ashape that corresponds to the refractive treatment shape. Further, thetreatment shape may be applied by providing the patient with a spectaclelens that has a shape that corresponds to the refractive treatmentshape. What is more, the treatment shape can be applied by providing thepatient with an intra-ocular lens that has a shape that corresponds tothe refractive treatment shape.

In another aspect, the present invention provides a system forgenerating a refractive treatment shape for ameliorating a visioncondition in an eye of a patient. The system comprises an input modulethat accepts a measurement surface aberration, a transformation modulethat derives a treatment surface aberration based on the measurementsurface aberration, and an output module that generates the refractivetreatment shape based on the treatment surface aberration. Themeasurement surface aberration may be based on a wavefront aberration ofthe eye. The wavefront aberration can correspond to a measurementsurface that is disposed at or near a pupil plane of the eye. Thetreatment surface aberration can correspond to a treatment surface thatis disposed at or near an anterior corneal surface of the eye, or atreatment surface that corresponds to a spectacle plane of the eye.Relatedly, the treatment surface may be disposed posterior to a pupilplane of the eye. The treatment surface aberration may be based on adifference between the measurement surface and the treatment surface.

In another aspect, the present invention provides a system forameliorating a vision condition in an eye of a patient. The systemcomprises an input module that accepts a measurement surface aberration,a transformation module that derives a treatment surface aberrationbased on the measurement surface aberration, an output module thatgenerates a refractive treatment shape based on the treatment surfaceaberration, and a laser system that directs laser energy onto the eyeaccording to the refractive treatment shape so as to reprofile a surfaceof the eye from an initial shape to a subsequent shape, the subsequentshape having correctively improved optical properties for amelioratingthe vision condition. The measurement surface aberration may be based ona wavefront aberration of the eye. The wavefront aberration cancorrespond to a measurement surface that is disposed at or near a pupilplane of the eye, and the treatment surface aberration can correspond toa treatment surface that is disposed at or near an anterior surface of acornea of the eye. The treatment surface aberration can be derived basedon a difference between the measurement surface and the treatmentsurface.

In some aspects, the treatment surface aberration may be a treatmentsurface wavefront map. In other aspects, the measurement surfaceaberration may be a measurement surface wavefront map. The treatmentsurface wavefront map may be derived at least in part by local slopescaling of the measurement surface wavefront map. In still otheraspects, the treatment surface wavefront map may be derived at least inpart by applying a scaling factor of 1/(1+Pd) to a slope of themeasurement surface wavefront map, where P represents a local curvatureof the measurement surface wavefront map and d represents a differencebetween the measurement surface and the treatment surface. In a relatedaspect, a difference between the measurement surface and a retinalsurface of the eye corresponds to a first vertex measure, and adifference between the treatment surface and the retinal surface of theeye corresponds to a second vertex measure. P may be based on a secondderivative of the measurement surface wavefront map. P may also be basedon a pupil radius of the eye.

In some aspects, the treatment surface wavefront map can be derivedaccording to an iterative Fourier reconstruction algorithm. What ismore, the measurement surface aberration may reflect low order and/orhigh order aberrations of the eye of the patient.

In another aspect, the present invention provides a system forgenerating a prescription for ameliorating a vision condition in an eyeof a patient. The system comprises an input that accepts irregularaberration data corresponding to an aberration measurement surfaceadjacent a pupil plane of the eye, a transformation module that derivesa treatment surface aberration corresponding to a treatment surface thatis disposed adjacent an anterior surface of a cornea of the eye, and anoutput module that generates the prescription based on the treatmentsurface aberration. The treatment surface aberration can be derived fromthe irregular aberration data using a difference between the measurementsurface and the treatment surface.

In one aspect, embodiments of the present invention provide a method ofcalculating a refractive treatment shape for ameliorating a visioncondition in an eye of a patient. The method can include determining afirst wavefront measurement corresponding to the pupil plane of the eye.The first wavefront measurement can include a first wavefront boundaryand a first set of wavefront coefficients. The method can also includedetermining a propagation distance between the pupil plane of the eyeand a treatment surface, and determining a propagated wavefrontmeasurement corresponding to the treatment surface based on the firstwavefront measurement and the propagation distance. The propagatedwavefront measurement can include a second wavefront boundary and asecond set of wavefront coefficients. The method can also includecalculating the refractive treatment shape based on the propagatedwavefront measurement. In some cases, the treatment surface correspondsto a corneal surface, a spectacle surface, a scleral lens surface, acontact lens surface, or an intraocular lens surface.

In another aspect, embodiments of the present invention provide a methodof ameliorating a vision condition in an eye of a patient. The methodincludes determining a first wavefront measurement corresponding to thepupil plane of the eye. The first wavefront measurement can include afirst wavefront boundary and a first set of wavefront coefficients. Themethod can also include determining a propagation distance between thepupil plane of the eye and a treatment surface, and determining apropagated wavefront measurement corresponding to the treatment surfacebased on the first wavefront measurement and the propagation distance.The propagated wavefront measurement can include a second wavefrontboundary and a second set of wavefront coefficients. The method can alsoinclude calculating the refractive treatment shape based on thepropagated wavefront measurement, and applying the refractive treatmentshape to the eye of the patient to ameliorate the vision condition. Insome cases, the refractive treatment shape is applied to the eye of thepatient in a treatment modality such as ablating a corneal surface ofthe eye to provide a corneal surface shape that corresponds to therefractive treatment shape, providing the patient with a contact lensthat has a shape that corresponds to the refractive treatment shape,providing the patient with a spectacle that has a shape that correspondsto the refractive treatment shape, providing the patient with a sclerallens that has a shape that corresponds to the refractive treatmentshape, or providing the patient with an intraocular lens that has ashape that corresponds to the refractive treatment shape.

In still another aspect, embodiments of the present invention provide asystem for generating a refractive treatment shape for ameliorating avision condition in an eye of a patient. The system can include an inputmodule that accepts a first wavefront measurement corresponding to thepupil plane of the eye. The first wavefront measurement can include afirst wavefront boundary and a first set of wavefront coefficients. Thesystem can also include a transformation module that derives apropagated wavefront measurement. The propagated wavefront measurementcan correspond to a treatment surface of the eye and include a secondwavefront boundary and a second set of wavefront coefficients. Thepropagated wavefront measurement can be derived from the first wavefrontmeasurement using a propagation distance between the first wavefrontmeasurement and the treatment surface. The system can also include anoutput module that generates the refractive treatment shape based on thepropagated wavefront measurement.

These and other aspects will be apparent in the remainder of thefigures, description, and claims.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 illustrates a laser ablation system according to an embodiment ofthe present invention.

FIG. 2 illustrates a simplified computer system according to anembodiment of the present invention.

FIG. 3 illustrates a wavefront measurement system according to anembodiment of the present invention.

FIG. 3A illustrates another wavefront measurement system according to anembodiment of the present invention.

FIG. 4 schematically represents a simplified set of modules that carryout one method of the present invention.

FIG. 5 is a flow chart that schematically illustrates a method ofdetermining a refractive treatment shape according to one embodiment ofthe present invention.

FIG. 5A depicts aspects of ocular wavefront sensing for Hartmann-Shackaberrometry according to embodiments of the present invention.

FIG. 6 illustrates a model optical system according to an embodiment ofthe present invention.

FIG. 6A shows geometry for a vertex correction for myopic and hyperopiccases according to embodiments of the present invention.

FIGS. 7A, 7B and 7C illustrate a comparison between vertex correctedpower calculations based on algorithms provided by the present inventionwith calculations based on a classical formula according to anembodiment of the present invention.

FIG. 8 illustrates a wavefront before and after a vertex correctionaccording to an embodiment of the present invention.

FIG. 8A shows examples of a diverging wavefront and a wavefront with aspherical aberration before and after propagation according toembodiments of the present invention.

FIG. 8B shows a geometry of a myopic correction as an original wavefrontand a propagated wavefront according to embodiments of the presentinvention.

FIG. 8C shows a geometry of a hyperopic wavefront and a correspondingpropagated wavefront according to embodiments of the present invention.

FIG. 9 shows a wavefront coordinate system according to an embodiment ofthe present invention.

FIG. 10 shows a wavefront coordinate system according to an embodimentof the present invention.

FIG. 11 depicts a wavefront coordinate system according to an embodimentof the present invention.

FIG. 12 depicts a wavefront coordinate system according to an embodimentof the present invention.

FIG. 12A shows how a wavefront boundary can change when the wavefrontpropagates, according to embodiments of the present invention.

FIG. 13 illustrates a procedural flow chart for analyzing wavefrontpropagation according to an embodiment of the present invention.

FIG. 14 illustrates a procedural flow chart for analyzing wavefrontpropagation according to an embodiment of the present invention.

FIG. 15 shows a procedural flow chart for analyzing wavefrontpropagation according to an embodiment of the present invention.

DETAILED DESCRIPTION OF THE INVENTION

The present invention provides methods, software, and systems forprocessing optical power measurements taken at a first plane andconverting those power measurements to corrective optical powermeasurements that can be used at a second plane.

The present invention is generally useful for enhancing the accuracy andefficacy of laser eye surgical procedures, such as photorefractivekeratectomy (PRK), phototherapeutic keratectomy (PTK), laser in situkeratomileusis (LASIK), and the like. The present invention can provideenhanced optical accuracy of refractive procedures by improving themethodology for processing measured optical errors of the eye and hencecalculate a more accurate refractive ablation program. In one particularembodiment, the present invention is related to therapeuticwavefront-based ablations of pathological eyes.

The present invention can be readily adapted for use with existing lasersystems, wavefront measurement systems, and other optical measurementdevices. While the systems, software, and methods of the presentinvention are described primarily in the context of a laser eye surgerysystem, it should be understood the present invention may be adapted foruse in alternative eye treatment procedures, systems, or modalities,such as spectacle lenses, intraocular lenses, contact lenses, cornealring implants, collagenous corneal tissue thermal remodeling, cornealinlays, corneal onlays, other corneal implants or grafts, and the like.Relatedly, systems, software, and methods according to embodiments ofthe present invention are well suited for customizing any of thesetreatment modalities to a specific patient. Thus, for example,embodiments encompass custom intraocular lenses, custom contact lenses,custom corneal implants, and the like, which can be configured to treator ameliorate any of a variety of vision conditions in a particularpatient based on their unique ocular characteristics or anatomy.

Turning now to the drawings, FIG. 1 illustrates a laser eye surgerysystem 10 of the present invention, including a laser 12 that produces alaser beam 14. Laser 12 is optically coupled to laser delivery optics16, which directs laser beam 14 to an eye E of patient P. A deliveryoptics support structure (not shown here for clarity) extends from aframe 18 supporting laser 12. A microscope 20 is mounted on the deliveryoptics support structure, the microscope often being used to image acornea of eye E.

Laser 12 generally comprises an excimer laser, ideally comprising anargon-fluorine laser producing pulses of laser light having a wavelengthof approximately 193 nm. Laser 12 will preferably be designed to providea feedback stabilized fluence at the patient's eye, delivered viadelivery optics 16. The present invention may also be useful withalternative sources of ultraviolet or infrared radiation, particularlythose adapted to controllably ablate the corneal tissue without causingsignificant damage to adjacent and/or underlying tissues of the eye.Such sources include, but are not limited to, solid state lasers andother devices which can generate energy in the ultraviolet wavelengthbetween about 185 and 215 nm and/or those which utilizefrequency-multiplying techniques. Hence, although an excimer laser isthe illustrative source of an ablating beam, other lasers may be used inthe present invention.

Laser system 10 will generally include a computer or programmableprocessor 22. Processor 22 may comprise (or interface with) aconventional PC system including the standard user interface devicessuch as a keyboard, a display monitor, and the like. Processor 22 willtypically include an input device such as a magnetic or optical diskdrive, an internet connection, or the like. Such input devices willoften be used to download a computer executable code from a tangiblestorage media 29 embodying any of the methods of the present invention.Tangible storage media 29 may take the form of a floppy disk, an opticaldisk, a data tape, a volatile or non-volatile memory, RAM, or the like,and the processor 22 will include the memory boards and other standardcomponents of modern computer systems for storing and executing thiscode. Tangible storage media 29 may optionally embody wavefront sensordata, wavefront gradients, a wavefront elevation map, a treatment map, acorneal elevation map, and/or an ablation table. While tangible storagemedia 29 will often be used directly in cooperation with a input deviceof processor 22, the storage media may also be remotely operativelycoupled with processor by means of network connections such as theinternet, and by wireless methods such as infrared, Bluetooth, or thelike.

Laser 12 and delivery optics 16 will generally direct laser beam 14 tothe eye of patient P under the direction of a computer 22. Computer 22will often selectively adjust laser beam 14 to expose portions of thecornea to the pulses of laser energy so as to effect a predeterminedsculpting of the cornea and alter the refractive characteristics of theeye. In many embodiments, both laser beam 14 and the laser deliveryoptical system 16 will be under computer control of processor 22 toeffect the desired laser sculpting process, with the processor effecting(and optionally modifying) the pattern of laser pulses. The pattern ofpulses may by summarized in machine readable data of tangible storagemedia 29 in the form of a treatment table, and the treatment table maybe adjusted according to feedback input into processor 22 from anautomated image analysis system in response to feedback data providedfrom an ablation monitoring system feedback system. Optionally, thefeedback may be manually entered into the processor by a systemoperator. Such feedback might be provided by integrating the wavefrontmeasurement system described below with the laser treatment system 10,and processor 22 may continue and/or terminate a sculpting treatment inresponse to the feedback, and may optionally also modify the plannedsculpting based at least in part on the feedback. Measurement systemsare further described in U.S. Pat. No. 6,315,413, the full disclosure ofwhich is incorporated herein by reference.

Laser beam 14 may be adjusted to produce the desired sculpting using avariety of alternative mechanisms. The laser beam 14 may be selectivelylimited using one or more variable apertures. An exemplary variableaperture system having a variable iris and a variable width slit isdescribed in U.S. Pat. No. 5,713,892, the full disclosure of which isincorporated herein by reference. The laser beam may also be tailored byvarying the size and offset of the laser spot from an axis of the eye,as described in U.S. Pat. Nos. 5,683,379, 6,203,539, and 6,331,177, thefull disclosures of which are incorporated herein by reference.

Still further alternatives are possible, including scanning of the laserbeam over the surface of the eye and controlling the number of pulsesand/or dwell time at each location, as described, for example, by U.S.Pat. No. 4,665,913, the full disclosure of which is incorporated hereinby reference; using masks in the optical path of laser beam 14 whichablate to vary the profile of the beam incident on the cornea, asdescribed in U.S. Pat. No. 5,807,379, the full disclosure of which isincorporated herein by reference; hybrid profile-scanning systems inwhich a variable size beam (typically controlled by a variable widthslit and/or variable diameter iris diaphragm) is scanned across thecornea; or the like. The computer programs and control methodology forthese laser pattern tailoring techniques are well described in thepatent literature.

Additional components and subsystems may be included with laser system10, as should be understood by those of skill in the art. For example,spatial and/or temporal integrators may be included to control thedistribution of energy within the laser beam, as described in U.S. Pat.No. 5,646,791, the full disclosure of which is incorporated herein byreference. Ablation effluent evacuators/filters, aspirators, and otherancillary components of the laser surgery system are known in the art.Further details of suitable systems for performing a laser ablationprocedure can be found in commonly assigned U.S. Pat. Nos. 4,665,913,4,669,466, 4,732,148, 4,770,172, 4,773,414, 5,207,668, 5,108,388,5,219,343, 5,646,791 and 5,163,934, the complete disclosures of whichare incorporated herein by reference. Basis data can be furthercharacterized for particular lasers or operating conditions, by takinginto account localized environmental variables such as temperature,humidity, airflow, and aspiration.

FIG. 2 is a simplified block diagram of an exemplary computer system 22that may be used by the laser surgical system 10 of the presentinvention. Computer system 22 typically includes at least one processor52 which may communicate with a number of peripheral devices via a bussubsystem 54. These peripheral devices may include a storage subsystem56, comprising a memory subsystem 58 and a file storage subsystem 60,user interface input devices 62, user interface output devices 64, and anetwork interface subsystem 66. Network interface subsystem 66 providesan interface to outside networks 68 and/or other devices, such as thewavefront measurement system 30.

User interface input devices 62 may include a keyboard, pointing devicessuch as a mouse, trackball, touch pad, or graphics tablet, a scanner,foot pedals, a joystick, a touchscreen incorporated into the display,audio input devices such as voice recognition systems, microphones, andother types of input devices. User input devices 62 will often be usedto download a computer executable code from a tangible storage media 29embodying any of the methods of the present invention. In general, useof the term “input device” is intended to include a variety ofconventional and proprietary devices and ways to input information intocomputer system 22.

User interface output devices 64 may include a display subsystem, aprinter, a fax machine, or non-visual displays such as audio outputdevices. The display subsystem may be a cathode ray tube (CRT), aflat-panel device such as a liquid crystal display (LCD), a projectiondevice, or the like. The display subsystem may also provide a non-visualdisplay such as via audio output devices. In general, use of the term“output device” is intended to include a variety of conventional andproprietary devices and ways to output information from computer system22 to a user.

Storage subsystem 56 can store the basic programming and data constructsthat provide the functionality of the various embodiments of the presentinvention. For example, a database and modules implementing thefunctionality of the methods of the present invention, as describedherein, may be stored in storage subsystem 56. These software modulesare generally executed by processor 52. In a distributed environment,the software modules may be stored on a plurality of computer systemsand executed by processors of the plurality of computer systems. Storagesubsystem 56 typically comprises memory subsystem 58 and file storagesubsystem 60.

Memory subsystem 58 typically includes a number of memories including amain random access memory (RAM) 70 for storage of instructions and dataduring program execution and a read only memory (ROM) 72 in which fixedinstructions are stored. File storage subsystem 60 provides persistent(non-volatile) storage for program and data files, and may includetangible storage media 29 (FIG. 1) which may optionally embody wavefrontsensor data, wavefront gradients, a wavefront elevation map, a treatmentmap, and/or an ablation table. File storage subsystem 60 may include ahard disk drive, a floppy disk drive along with associated removablemedia, a Compact Digital Read Only Memory (CD-ROM) drive, an opticaldrive, DVD, CD-R, CD-RW, solid-state removable memory, and/or otherremovable media cartridges or disks. One or more of the drives may belocated at remote locations on other connected computers at other sitescoupled to computer system 22. The modules implementing thefunctionality of the present invention may be stored by file storagesubsystem 60.

Bus subsystem 54 provides a mechanism for letting the various componentsand subsystems of computer system 22 communicate with each other asintended. The various subsystems and components of computer system 22need not be at the same physical location but may be distributed atvarious locations within a distributed network. Although bus subsystem54 is shown schematically as a single bus, alternate embodiments of thebus subsystem may utilize multiple busses.

Computer system 22 itself can be of varying types including a personalcomputer, a portable computer, a workstation, a computer terminal, anetwork computer, a control system in a wavefront measurement system orlaser surgical system, a mainframe, or any other data processing system.Due to the ever-changing nature of computers and networks, thedescription of computer system 22 depicted in FIG. 2 is intended only asa specific example for purposes of illustrating one embodiment of thepresent invention. Many other configurations of computer system 22 arepossible having more or less components than the computer systemdepicted in FIG. 2.

Referring now to FIG. 3, one embodiment of a wavefront measurementsystem 30 is schematically illustrated in simplified form. In verygeneral terms, wavefront measurement system 30 is configured to senselocal slopes of a gradient map exiting the patient's eye. Devices basedon the Hartmann-Shack principle generally include a lenslet array tosample the gradient map uniformly over an aperture, which is typicallythe exit pupil of the eye. Thereafter, the local slopes of the gradientmap are analyzed so as to reconstruct the wavefront surface or map.

More specifically, one wavefront measurement system 30 includes an imagesource 32, such as a laser, which projects a source image throughoptical tissues 34 of eye E so as to form an image 44 upon a surface ofretina R. The image from retina R is transmitted by the optical systemof the eye (e.g., optical tissues 34) and imaged onto a wavefront sensor36 by system optics 37. The wavefront sensor 36 communicates signals toa computer system 22′ for measurement of the optical errors in theoptical tissues 34 and/or determination of an optical tissue ablationtreatment program. Computer 22′ may include the same or similar hardwareas the computer system 22 illustrated in FIGS. 1 and 2. Computer system22′ may be in communication with computer system 22 that directs thelaser surgery system 10, or some or all of the components of computersystem 22, 22′ of the wavefront measurement system 30 and laser surgerysystem 10 may be combined or separate. If desired, data from wavefrontsensor 36 may be transmitted to a laser computer system 22 via tangiblemedia 29, via an I/O port, via an networking connection 66 such as anintranet or the Internet, or the like.

Wavefront sensor 36 generally comprises a lenslet array 38 and an imagesensor 40. As the image from retina R is transmitted through opticaltissues 34 and imaged onto a surface of image sensor 40 and an image ofthe eye pupil P is similarly imaged onto a surface of lenslet array 38,the lenslet array separates the transmitted image into an array ofbeamlets 42, and (in combination with other optical components of thesystem) images the separated beamlets on the surface of sensor 40.Sensor 40 typically comprises a charged couple device or “CCD,” andsenses the characteristics of these individual beamlets, which can beused to determine the characteristics of an associated region of opticaltissues 34. In particular, where image 44 comprises a point or smallspot of light, a location of the transmitted spot as imaged by a beamletcan directly indicate a local gradient of the associated region ofoptical tissue.

Eye E generally defines an anterior orientation ANT and a posteriororientation POS. Image source 32 generally projects an image in aposterior orientation through optical tissues 34 onto retina R asindicated in FIG. 3. Optical tissues 34 again transmit image 44 from theretina anteriorly toward wavefront sensor 36. Image 44 actually formedon retina R may be distorted by any imperfections in the eye's opticalsystem when the image source is originally transmitted by opticaltissues 34. Optionally, image source projection optics 46 may beconfigured or adapted to decrease any distortion of image 44.

In some embodiments, image source optics 46 may decrease lower orderoptical errors by compensating for spherical and/or cylindrical errorsof optical tissues 34. Higher order optical errors of the opticaltissues may also be compensated through the use of an adaptive opticelement, such as a deformable mirror (described below). Use of an imagesource 32 selected to define a point or small spot at image 44 uponretina R may facilitate the analysis of the data provided by wavefrontsensor 36. Distortion of image 44 may be limited by transmitting asource image through a central region 48 of optical tissues 34 which issmaller than a pupil 50, as the central portion of the pupil may be lessprone to optical errors than the peripheral portion. Regardless of theparticular image source structure, it will be generally be beneficial tohave a well-defined and accurately formed image 44 on retina R.

In one embodiment, the wavefront data may be stored in a computerreadable medium 29 or a memory of the wavefront sensor system 30 in twoseparate arrays containing the x and y wavefront gradient valuesobtained from image spot analysis of the Hartmann-Shack sensor images,plus the x and y pupil center of fsets from the nominal center of theHartmann-Shack lenslet array, as measured by the pupil camera 51 (FIG.3) image. Such information contains all the available information on thewavefront error of the eye and is sufficient to reconstruct thewavefront or any portion of it. In such embodiments, there is no need toreprocess the Hartmann-Shack image more than once, and the data spacerequired to store the gradient array is not large. For example, toaccommodate an image of a pupil with an 8 mm diameter, an array of a20×20 size (i.e., 400 elements) is often sufficient. As can beappreciated, in other embodiments, the wavefront data may be stored in amemory of the wavefront sensor system in a single array or multiplearrays.

While the methods of the present invention will generally be describedwith reference to sensing of an image 44, it should be understood that aseries of wavefront sensor data readings may be taken. For example, atime series of wavefront data readings may help to provide a moreaccurate overall determination of the ocular tissue aberrations. As theocular tissues can vary in shape over a brief period of time, aplurality of temporally separated wavefront sensor measurements canavoid relying on a single snapshot of the optical characteristics as thebasis for a refractive correcting procedure. Still further alternativesare also available, including taking wavefront sensor data of the eyewith the eye in differing configurations, positions, and/ororientations. For example, a patient will often help maintain alignmentof the eye with wavefront measurement system 30 by focusing on afixation target, as described in U.S. Pat. No. 6,004,313, the fulldisclosure of which is incorporated herein by reference. By varying aposition of the fixation target as described in that reference, opticalcharacteristics of the eye may be determined while the eye accommodatesor adapts to image a field of view at a varying distance and/or angles.

The location of the optical axis of the eye may be verified by referenceto the data provided from a pupil camera 52. In the exemplaryembodiment, a pupil camera 52 images pupil 50 so as to determine aposition of the pupil for registration of the wavefront sensor datarelative to the optical tissues.

An alternative embodiment of a wavefront measurement system isillustrated in FIG. 3A. The major components of the system of FIG. 3Aare similar to those of FIG. 3. Additionally, FIG. 3A includes anadaptive optical element 53 in the form of a deformable mirror. Thesource image is reflected from deformable mirror 98 during transmissionto retina R, and the deformable mirror is also along the optical pathused to form the transmitted image between retina R and imaging sensor40. Deformable mirror 98 can be controllably deformed by computer system22 to limit distortion of the image formed on the retina or ofsubsequent images formed of the images formed on the retina, and mayenhance the accuracy of the resultant wavefront data. The structure anduse of the system of FIG. 3A are more fully described in U.S. Pat. No.6,095,651, the full disclosure of which is incorporated herein byreference.

The components of an embodiment of a wavefront measurement system formeasuring the eye and ablations comprise elements of a VISX WaveScan®,available from VISX, INCORPORATED of Santa Clara, Calif. One embodimentincludes a WaveScan® with a deformable mirror as described above. Analternate embodiment of a wavefront measuring system is described inU.S. Pat. No. 6,271,915, the full disclosure of which is incorporatedherein by reference.

FIG. 4 schematically illustrates a simplified set of modules, or acorrection system 100, for carrying out a method according to oneembodiment of the present invention. Correction system 100 can beintegrated or interfaced with, for example, computer system 22, orotherwise used in conjunction with laser surgical system 10. The modulesmay be software modules on a computer readable medium that is processedby processor 52 (FIG. 2), hardware modules, or a combination thereof.Any of a variety of commonly used platforms, such as Windows, MacIntosh,and Unix, along with any of a variety of commonly used programminglanguages, may be used to implement the present invention.

Correction system 100 can be configured to generate a refractivetreatment shape 110 for ameliorating a vision condition in a patient. Aninput module 102 typically receives a measurement surface aberration120, such as wavefront aberration data from wavefront sensors, whichcharacterize aberrations and other optical characteristics of the entireoptical tissue system imaged. Often, the wavefront aberrationcorresponds to a measurement surface that is disposed at or near a pupilplane of the eye. The data from the wavefront sensors are typicallygenerated by transmitting an image (such as a small spot or point oflight) through the optical tissues, as described above. Measurementsurface aberration 120 can include an array of optical gradients or agradient map.

Correction system 100 can include a transformation module 104 thatderives a treatment surface aberration. The treatment surface aberrationcan correspond to a treatment surface that is disposed at or near ananterior corneal surface of the eye, or a treatment surface thatcorresponds to a spectacle plane of the eye. Relatedly, the treatmentsurface may be disposed posterior to a pupil plane of the eye. Often,the treatment surface aberration is derived from measurement surfaceaberration 120 using a difference between the measurement surface andthe treatment surface. For example, optical gradient data from inputmodule 102 may be transmitted to transformation module 104, where atreatment surface aberration is mathematically reconstructed based onthe optical gradient data.

Correction system 100 can include an output module 106, such that thetreatment surface aberration generated by transformation module 104 canthen be transmitted to output module 106 where a refractive treatmentshape 110 can be generated based on the treatment surface aberration.Refractive treatment shape 110 may be transmitted to a laser treatmentapparatus for generation of a laser ablation treatment for the patient.Similarly, refractive treatment shape 110 may form the basis forfabrication of contact lenses, spectacles, or intra-ocular lenses.

As can be appreciated, the present invention should not be limited tothe order of steps, or the specific steps illustrated, and variousmodifications to the method, such as having more or less steps, may bemade without departing from the scope of the present invention.

In one embodiment, the present invention provides a method ofdetermining a refractive treatment shape for ameliorating a visioncondition in a patient. FIG. 5 depicts the steps of an exemplary methodaccording to the present invention. The refractive treatment shape canbe based on a treatment surface aberration that is derived from ameasurement surface aberration.

For ocular wavefront measurements, the Hartmann-Shack wavefront sensor(Liang, J. et al., J. Opt. Soc. Am. A 11:1949-1957 (1994)) can be usedas an aberrometer. FIGS. 5A(a) and 5A(b) show aspects of the ocularwavefront sensing for the Hartmann-Shack aberrometry according toembodiments of the present invention. FIG. 5A(a) shows that a parallelwave can be formed when the rays of a beacon from the retina passthrough an eye with no ocular aberrations. As shown in FIG. 5A(b), whenthe optics of the entire eye are simplified as a virtual lens, theeffective pupil can be one that is magnified and anterior to the virtuallens. As shown in FIG. 5A(a), a thin beam of light can be projected ontothe retina to form a beacon in such an outgoing wavefront sensingsystem. The rays of the pseudo point source pass through the optics ofthe eye to become a plane wave, if there are no ocular aberrations inthe eye. However, if ocular aberrations exist, the rays would not form aplane wave and the deviations in the optical path difference can bemeasured as the ocular wavefront.

The ocular wavefront diameter can be determined by the pupil size.Because of the large dioptric power of the cornea, a significantmagnification of the pupil size may occur. In the optics of aHartmann-Shack wavefront sensor, the detected pupil size on a CCD camerais typically the magnified pupil size, not the physical pupil size. Therepresentation of a wavefront may therefore use a magnified pupil size,as can be explained by FIG. 5A(b). In FIG. 5A(b), the entire optics ofthe eye are simplified as a thin lens, and the effective pupil isimmediately anterior to the simplified lens. Because the parallel waveis formed after it leaves the entire optics of the eye, the effectivepupil size is tied to the size when the wavefront becomes parallel.Therefore, the magnified pupil size can be taken as the wavefrontdiameter.

However, the ocular wavefront is often measured on the exit pupil plane.Even though the boundary of the wavefront is the magnified pupil size,it is different than when it propagates to the curved cornea surface orother planes, such as a spectacle plane. The problem for the powerchange for correcting astigmatic eyes at different spectacle planes wasinvestigated extensively by Harris (Harris, W. F. Optom. Vis. Sci.73:606-612 (1996)). For ocular wavefronts that include high orderaberrations, the problem for such propagations has not been considered.Although the light rays form a diverging wave before they leave thesurface of the cornea, they become parallel wave when they pass thesurface of the cornea. Embodiments of the present invention address freespace wavefront propagation with the flat surface as a reference plane.Embodiments also encompass the propagation of the spherocylindricalocular aberrations, as well as a more generic approach that includes thetreatment of the wavefront propagation of both low order and high orderaberrations.

I. Measurement Surface Aberration

In general terms, a measurement surface aberration can be determinedfrom optical data corresponding to a measurement surface. For example, ameasurement surface aberration can be determined by measuring awavefront aberration of an eye of a patient. Measurement surfaceaberrations can be determined by aberrometers such as Hartmann-Shackaberrometers, ray tracing aberrometers, Tscherning aberrometers,Scheiner aberrometers, double-pass aberrometers, and the like, as wellas topographical devices. In some embodiments, a wavefront measurementsystem that includes a wavefront sensor (such as a Hartmann-Shacksensor) may be used to obtain one or more measurement surfaceaberrations (e.g. wavefront maps) of the optical tissues of the eye. Thewavefront map may be obtained by transmitting an image through theoptical tissues of the eye and sensing the exiting wavefront surface.From the wavefront map, it is possible to calculate a surface gradientor gradient map across the optical tissues of the eye. A gradient mapmay comprise an array of the localized gradients as calculated from eachlocation for each lenslet of the Hartmann-Shack sensor. Measurementsurface aberrations can correspond to a first plane or surface, and canencompass a first measurement surface boundary and a first measurementsurface set of coefficients. For example, a first wavefront measurementcan include a first wavefront boundary and a first set of wavefrontcoefficients. In some embodiments, a system or method may involve ameasurement surface aberration corresponding to a measurement surface ofthe eye. The measurement surface aberration may include a measurementsurface boundary and a measurement surface magnitude.

A. Measurement Surface

There are a variety of devices and methods for measuring surfacecharacteristics of optical systems. The category of aberroscopes oraberrometers includes classical phoropter and wavefront approaches.Classical phoropters can be used to record optical data corresponding toa measurement surface that is disposed anterior to the cornea of an eye.For example, phoropters can measure low order aberrations at a distanceof about 12.5 mm anterior to the corneal surface. In many cases, thiswill correspond to a spectacle plane of the eye. Wavefront devices areoften used to measure both low order and high order aberrations adjacenta pupil plane, which can be about 3.5 mm posterior to the cornealsurface. Another category of measuring approaches includes topographybased measuring devices and methods. Topography typically involvesproviding optical data corresponding to a measurement surface that isdisposed at or near the corneal surface of the eye. In some embodiments,the terms “plane” and “surface” may be used interchangeably.

B. Aberration

As noted above, the measurement surface aberration can be based on arefractive measurement as determined by an optometer, or any of a widevariety of devices for obtaining irregular aberration data. Similarly,the measurement surface aberration can be a measurement surfacewavefront map, as determined by a wavefront measurement device. What ismore, the measurement surface aberration may reflect both low order andhigh order aberrations of the eye of a patient. In some cases,aberrations can be embodied as an elevation map, a surface map, or thelike.

II. Treatment Surface Aberration

When a measurement surface aberration of an optical system has beendetermined, it is then possible to derive a treatment surface aberrationof the optical system. In the case of refractive surgical methods, forexample, a treatment surface aberration corresponding to a corneal planecan be derived from a measurement surface aberration as determined in aplane other than the corneal plane. Treatment surface aberrations cancorrespond to a second plane or surface, and can encompass a secondsurface boundary and a second surface set of coefficients. For example,a treatment surface aberration can include a first treatment boundaryand a first set of treatment coefficients. In some embodiments, systemsor methods can involve a propagation distance between a measurementsurface of the eye and a treatment surface of the eye. Embodiments alsoencompass a treatment surface aberration based on the measurementsurface aberration and the propagation distance. A treatment surfaceaberration can include a treatment surface boundary and a treatmentsurface magnitude.

A. Treatment Surface

The treatment surface aberration corresponds to a treatment surface,which is typically disposed at or near an anterior surface of a corneaof an eye. Often, the treatment surface will correspond to a cornealplane associated with the eye, as in the case of ablative laser eyesurgery or contact lens treatments. At other times, the treatmentsurface may correspond to a spectacle plane associated with the eye, asin the case of spectacle treatments. Further, the treatment surface canbe posterior to the pupil plane of the eye, as in the case ofintraocular lens treatments. As noted above, in some embodiments, theterms “plane” and “surface” may be used interchangeably.

B. Derivation of Treatment Surface Aberration

The treatment surface aberration can be derived from the measurementsurface aberration, based on a difference between the measurementsurface and the treatment surface. The difference between themeasurement surface and the treatment surface, for example, can includea distance measurement that represents a difference between the twosurfaces. In some embodiments, the distance measurement is based on avertex distance difference, the vertex distance difference reflecting adistance between a vertex of the measurement surface and a vertex of thetreatment surface.

1. Classical Vertex Correction Formulas

Traditionally, the power of a lens is measured in diopters, and can bedefined as the reciprocal of the lens focal length in meters. FIG. 6shows a schematic diagram of an optical system. The system includes afirst plane disposed at a first distance from a focal plane, the firstdistance corresponding to a first focal length, and a second planedisposed at a second distance from the focal plane, the second distancecorresponding to a second focal length. Although the first and secondplanes are illustrated as flat surfaces, these planes can also representcurved surfaces such as lenses, wavefronts, and other representations ofoptical surfaces or systems. In the exemplary optical system depicted byFIG. 6 legend (a), the focal plane may be associated with a retinalplane, the first plane may be associated with a spectacle plane, and thesecond plane may be associated with a corneal plane.

A treatment surface can correspond to, or be based upon, a spectaclesurface, corneal surface, pupil surface, and the like. A spectaclesurface is typically about 12.5 mm anterior to the cornea of the eye. Apupil surface or plane is typically about 3.5 mm posterior to the corneaof the eye. An intraocular lens surface is usually about 0.5 mmposterior to the pupil surface or plane of the eye. A contact lenssurface is typically at or slightly anterior to the cornea of the eye.If the treatment surface and the measurement surface are substantiallyin the same plane, there may be no need for a vertex correction.

When prescribing spectacles, for example, an optometrist may first makeor consider an aberration measurement such as a refractive measurementof the eye, where the aberration measurement corresponds to ameasurement surface at or near a pupil plane or surface of the eye.Because the treatment surface may not be the same as the measurementsurface, it is often desirable to make a power adjustment in order todetermine the corrective surface shape or treatment surface aberration.In the case of spectacles, the treatment surface is disposed anterior tothe corneal surface, usually by about 12.5 mm.

Likewise, when prescribing contact lenses, an optometrist can consider arefractive correction corresponding to the spectacle plane, and make apower adjustment to account for the difference between the spectacleplane and the corneal plane. In this case, the adjustment often dependson a vertex distance, corresponding to the distance between theposterior surface of the spectacle lens and the anterior surface of thecornea.

Thus, given a measurement surface aberration, it is possible to derive atreatment surface aberration based on a difference between the treatmentsurface and the measurement surface. Often, the difference will be avertex distance between the treatment surface and the measurementsurface. As further discussed below, the treatment surface aberrationcan then be used to determine a refractive treatment shape. In the caseof corrective spectacles, the refractive treatment shape can be a basisfor a prescription for the patient, where the treatment shapecorresponds to the spectacle plane or surface.

Typically, the measurement surface aberration corresponds to a firstpower data, and the treatment surface aberration corresponds to aderived second power data. The second power data can be derived from thefirst power data and the distance between the measurement surface andthe treatment. To achieve the effect of a power change, in terms of avertex correction, a vertex distance measure can be based on adifference between the measurement surface and the treatment surface.The vertex correction represents a power change between the first powerdata and the second power data. In this sense, the derivation of thesecond power corresponds to a vertex correction of the first power. Thevertex of a lens curve can be defined as the apex of the lens curve, oras the point on the lens curve in which the lens curve axis intersectsit.

a. Traditional (Non Wavefront)

Traditional phoropters can be used to make traditional opticalaberration measurements such as sphere and cylinder. Such aberrationmeasurements are often expressed in terms of dioptric power. Referringagain to FIG. 6 legend (a), assuming the power corresponding to thesecond plane, e.g. a corneal plane, is S, and the power correspondingthe first plane, e.g. a spectacle plane, is S′, it is possible todescribe the relationship between the two powers with the followingequations.

$\begin{matrix}{{S = \frac{1}{f}},} & (1) \\{S^{\prime} = {\frac{1}{f + d} = \frac{S}{1 + {dS}}}} & (2)\end{matrix}$

Power can be expressed in units of diopters. f represents the distancebetween the focal plane and the second plane, although here this term isnot a factor in the relationship between the two power measurements Sand S′. d represents the vertex distance between the first and secondplanes. Where a first plane treatment surface is disposed anterior to asecond plane measurement surface, d will typically have a positivevalue. For example, for spectacle treatments, d can be about 0.0125 m,and for refractive surgery treatments, d can be about 0.0035 m.Conversely, where the first plane treatment shape is disposed posteriorto a second plane measurement surface, d will typically have a negativevalue. For example, for intraocular lens treatments, d can be about−0.0005 m.

Sphere is a low order aberration corresponding to defocus, and cylinderis a low order aberration corresponding to astigmatism. To consider acombination of sphere and cylinder powers, it is possible to replace Sby (S+C) where C stands for cylinder power at the maximum meridian.Thus, cylinder at the spectacle plane can be represented by C′, where

$\begin{matrix}{C^{\prime} = {\frac{S + C}{1 + {d\left( {S + C} \right)}} - {S^{\prime}.}}} & (3)\end{matrix}$

These formulae can be used to calculate the power change associated witha vertex distance.

A related embodiment that encompasses a formulation of the classicalvertex correction also addresses the situation where only low orderaberrations, namely, the sphere and cylinder, exist. FIG. 6A showsgeometries for the vertex correction for myopic and hyperopic casesaccording to embodiments of the present invention. FIG. 6A(a) shows ageometry for the vertex correction in a myopic case. FIG. 6A(b) shows ageometry for the vertex correction in a hyperopic case. In someembodiments, it is helpful to assume the original plane is less anteriorthan the new plane. In some cases, the focal length is fbefore thevertex correction and the vertex distance is d, both in meters. SupposeS stands for the sphere power and C stands for the cylinder power, bothin diopters, before the vertex correction. After the vertex correction,the sphere and cylinder are denoted as S′ and C′. Let's first considerthe pure sphere case. From the geometrical optics, for the myopic casewe have

$\begin{matrix}{{S = \frac{1}{f}},} & ({A1a}) \\{{S^{\prime} = \frac{S}{f - d^{\prime}}},} & ({A1b})\end{matrix}$

where f stands for the focal length and d for the vertex distance, bothin meters. Solving for f from Eq. (Ala) and substituting it into Eq.(A1b), we obtain the vertex correction for myopia as

$\begin{matrix}{S^{\prime} = {\frac{S}{1 - {Sd}}.}} & ({A2})\end{matrix}$

Similarly, for the hyperopic case, we have

$\begin{matrix}{{S = \frac{1}{f}},} & ({A3a}) \\{{S^{\prime} = \frac{1}{f + d}},} & ({A3b})\end{matrix}$

and the vertex correction formula for hyperopic can be obtained as

$\begin{matrix}{S^{\prime} = {\frac{S}{1 + {Sd}}.}} & ({A4})\end{matrix}$

Noticing that for myopia, the convention is that the sphere power isnegative so Eq. (A4) can be used for both cases of myopia and hyperopia.When the correction is in the reverse direction, i.e., changing thepower from a more anterior plane to a less anterior plane, the vertexdistance d should take a negative value.

For cylinder case, we only need to consider two meridians, the maximumpower and the minimum power before and after the vertex correction. Inplus cylinder notation, the maximum power is S+C and the minimum poweris S. In minus cylinder notation, the maximum power is S and the minimumpower is S+C. Therefore, only the two powers, S+C and S, need to bevertex corrected. With a similar approach, we can obtain the vertexcorrection formula for sphere and cylinder as

$\begin{matrix}{{S^{\prime} = \frac{S}{1 + {Sd}}},} & ({A5a}) \\{{S^{\prime} + C^{\prime}} = {\frac{S + C}{1 + {\left( {S + C} \right)d}}.}} & ({A5b})\end{matrix}$

Equation (A5) is the standard formula for vertex correction for loworder spherocylindrical error. Again, when the correction is changedfrom a more anterior plane to a less anterior plane, the vertex distanced should take a negative value.

b. Wavefront

In addition to the traditional phoropter approaches discussed above, itis also possible to evaluate optical systems based on wavefrontanalysis. Wavefront analysis can be useful in evaluating low order andhigh order aberrations. Referring again to FIG. 6, it is possible toconsider the first and second planes as associated with a generalwavefront. The wavefront can begin at a virtual focal pointcorresponding to the focal plane, and propagate from plane two towardplane one. For each point along the wavefront surface, a local slope canbe calculated. For example, the local slope can be the first derivativeof the surface at a certain point. The local slope reflects a surfacevalue at that point, as well as the surface values of the surroundingpoints. The local slope can be a direction-dependent vector. Because thewavefront local slopes are proportional to the local focal length, asthe wavefront is propagated forward, the slope of the wavefront can bescaled by a factor of α such that:

$\begin{matrix}{\alpha = \frac{f}{f + d}} & (4)\end{matrix}$

where f is the focal length of the wavefront and d is the vertexdistance. Here, the vertex distance can represent a difference betweenthe measurement surface, or plane two, and the treatment surface, orplane one. Thus, by making an initial measurement of the wavefront atplane two, it is possible to calculate a new wavefront surface at planewhere individual points on the new surface have a local curvature thatis derived by the scaling factor as discussed above. In the exemplaryoptical system depicted by FIG. 6 legend (b), the first plane canrepresent a corneal plane, the second plane can represent a pupil plane,and the focal plane can represent a retinal plane. If the treatmentsurface is anterior to the measurement surface, then the vertex distanceis positive, and if the treatment surface is posterior to themeasurement surface, then the vertex distance is negative. Similarly,for the myopia case, because the power is negative, the focal lengthcould take a negative value. Generally a can have a positive value, asthe absolute value of f is often much larger than d.

As discussed above, vertex correction can be used with traditionalaberrometry approaches. It is also possible to use vertex correctionwith wavefront approaches. Here, W(x,y) represents the wavefront at themeasurement plane and W′(x,y) represents the wavefront at the treatmentplane with vertex distance of d. The local slope is assumed to bescaled, as discussed above. Thus, the following equations are partialderivatives of the corrected wavefront at the treatment plane.

$\begin{matrix}{{\frac{\partial W^{\prime}}{\partial x} = {\frac{f}{f + d}\frac{\partial W}{\partial x}}}{\frac{\partial W^{\prime}}{\partial y} = {\frac{f}{f + d}\frac{\partial W}{\partial y}}}} & (5)\end{matrix}$

It can be demonstrated that the classical formula for vertex correctionholds with the assumption that the local slopes can be scaled accordingto a scaling factor of f/(f+d). The following examples illustrate thisprinciple with respect to (i) sphere, or defocus, (ii) cylinder, orastigmatism, (iii) coma, and (iv) spherical aberration. Wavefronts canbe expressed in terms of polynomial equations. This equation is usefulfor the derivations that follow.

$\begin{matrix}{\frac{\partial^{2}W^{\prime}}{\partial r^{2}} = {{\frac{x^{2}}{x^{2} + y^{2}}\frac{\partial^{2}W^{\prime}}{\partial x^{2}}} + {\frac{2{xy}}{x^{2} + y^{2}}\frac{\partial^{2}W^{\prime}}{{\partial x}{\partial y}}} + {\frac{y^{2}}{x^{2} + y^{2}}{\frac{\partial^{2}W^{\prime}}{\partial y^{2}}.}}}} & (6)\end{matrix}$

(i) Sphere

In the following discussion, Zernike polynomials are used to representthe ocular aberrations. Starting with a sphere, where W(r)=c₂ ⁰√{squareroot over (3)}(2r²−1), corresponding to the wavefront at the secondplane, the curvature of the converted wavefront W′(r) at the first planecan be expressed as

$\begin{matrix}{\begin{matrix}{\frac{\partial^{2}W^{\prime}}{\partial r^{2}} = {{\frac{x^{2}}{x^{2} + y^{2}}\frac{\partial^{2}W^{\prime}}{\partial x^{2}}} + {\frac{2{xy}}{x^{2} + y^{2}}\frac{\partial^{2}W^{\prime}}{{\partial x}{\partial y}}} + {\frac{y^{2}}{x^{2} + y^{2}}\frac{\partial^{2}W^{\prime}}{\partial y^{2}}}}} \\{{= {4\sqrt{3}c_{2}^{0}\frac{f}{f + d}}},}\end{matrix}{or}} & (7) \\{{\frac{\partial^{2}W^{\prime}}{\partial r^{2}} = {4\sqrt{3}c_{2}^{0}\frac{f}{f + d}}},} & (8)\end{matrix}$

where the curvature of the vertex corrected wavefront can be calculatedusing Equation (6). Here, r represents the normalized pupil radius withvalues from 0 to 1, x and y are the normalized values in x- and y-axis,f is the local focal length, or the reciprocal of local power, and c₂ ⁰is the Zernike coefficient of defocus term. From the definition ofpower, we have

$\begin{matrix}{{\frac{\partial^{2}W}{\partial r^{2}} = {4\sqrt{3}c_{2}^{0}}}{S = {\frac{1}{R^{2}}\frac{\partial^{2}W}{\partial r^{2}}}}{S^{\prime} = {\frac{1}{R^{2}}{\frac{\partial^{2}W^{\prime}}{\partial r^{2}}.}}}} & (9)\end{matrix}$

From Equations (8) and (9), we obtain the following formula

$\begin{matrix}{S^{\prime} = {{\frac{f}{f + d}S} = {\frac{S}{1 + {Sd}}.}}} & (10)\end{matrix}$

Equation (10) is the classical formula for vertex correction of puresphere power, thus demonstrating that vertex correction can beeffectively used in wavefront analysis.

(ii) Cylinder

In another example for astigmatism, W(r,θ)=c₂ ⁻²√{square root over(6)}r² sin 2θ+c₂ ²√{square root over (6)}r² cos 2θ corresponds to thewavefront at the second plane, a similar approach can be used to obtainthe curvature of the corrected wavefront as

$\begin{matrix}{\frac{\partial^{2}W^{\prime}}{\partial r^{2}} = {{\left( {{2\sqrt{6}c_{2}^{- 2}\sin \; 2\theta} + {2\sqrt{6}c_{2}^{2}\cos \; 2\theta}} \right)\frac{f}{f + d}} = {\frac{\partial^{2}W}{\partial r^{2}}{\frac{f}{f + d}.}}}} & (11)\end{matrix}$

Denoting P′ and P as the curvatures of W′ (converted wavefront) and W(measured wavefront) respectively,

$\begin{matrix}{P^{\prime} = {{P\frac{f}{f + d}} = {\frac{P}{1 + {Pd}}.}}} & (12)\end{matrix}$

By replacing P with S+C, it is possible to obtain the classical vertexcorrection for cylinder

$\begin{matrix}{C^{\prime} = {\frac{S + C}{1 + {d\left( {S + C} \right)}} - {S^{\prime}.}}} & (13)\end{matrix}$

(iii) Coma

In addition to the low order wavefront vertex corrections discussedabove, it is also possible to use vertex correction with wavefrontmeasurements that include high order aberrations. For example,horizontal coma can be expressed as W(r,θ)=√{square root over (8)}c₃¹(3r³−2r)cos θ, again corresponding to the wavefront at the secondplane. With an approach similar to that described above, it is possibleto calculate the derivatives to x and to y and then calculate thecurvature to r as

$\begin{matrix}{\frac{\partial^{2}W^{\prime}}{\partial r^{2}} = {{\frac{f}{f + d}18\sqrt{8}c_{3}^{1}x} = {\frac{\partial^{2}W}{\partial r^{2}}{\frac{f}{f + d}.}}}} & (14)\end{matrix}$

Denoting P′ and P as the curvatures of W′ (converted wavefront) and W(measured wavefront) respectively,

$\begin{matrix}{P^{\prime} = {{P\frac{f}{f + d}} = {\frac{P}{1 + {Pd}}.}}} & (15)\end{matrix}$

(iv) Spherical Aberrations

In another example, a spherical aberration can be expressed asW(r)=√{square root over (5)}c₄ ⁰(6r⁴−6r²+1). Again, with an approachsimilar to that described above, it is possible to calculate thederivatives to x and to y and then calculate the curvature to r todetermine the curvature of the corrected wavefront as

$\begin{matrix}{\frac{\partial^{2}W^{\prime}}{\partial r^{2}} = {{\frac{f}{f + d}\left( {{72r^{2}} - 12} \right)\sqrt{5}c_{4}^{0}} = {\frac{\partial^{2}W}{\partial r^{2}}{\frac{f}{f + d}.}}}} & (16)\end{matrix}$

Denoting P′ and P as the curvatures of W′ (converted wavefront) and W(measured wavefront) respectively,

$\begin{matrix}{P^{\prime} = {{P\frac{f}{f + d}} = {\frac{P}{1 + {Pd}}.}}} & (17)\end{matrix}$

Therefore, for low order aberrations as well as for high orderaberrations, it can be shown that by means of a slope scaling as appliedin wavefront, it is possible to achieve the effect of power change asdefined in a classical sense. Such approaches can be useful indetermining treatment surface aberrations based on measurement surfaceaberrations.

2. New Algorithm for Vertex Correction

Treatment surface aberrations can also be determined based on variousalgorithmic approaches. In some embodiments, the treatment surfaceaberration is a treatment surface wavefront map. The treatment surfacewavefront map can be derived at least in part by local slope scaling ofa measurement surface wavefront map. For example, a treatment surfacewavefront map can be derived at least in part by applying a scalingfactor of 1/(1+Pd) to a slope of a measurement surface wavefront map,where P represents a local curvature of the measurement surfacewavefront map and d represents a difference between a measurementsurface and a treatment surface. For example, P can be based on a secondderivative of the measurement surface wavefront map. P can also be basedon a pupil radius of the eye. The following examples illustratealgorithmic approaches that incorporate such principles.

a. Constant HOA

This algorithm assumes that the average curvature for low orderaberrations (LOA), as manifested by sphere and cylinder power terms, isaffected by vertex distance change. High order aberrations (HOA) areconsidered as local irregularity add-ons to the mean curvature, and arenot affected by vertex distance change. Thus, a total wavefront map canbe separated into low order and high order portions as shown by thefollowing formula

W(x,y)=W _(L)(x,y)+W _(H)(x,y).  (18)

For the low order portion, it is possible to obtain the sphere andcylinder power terms by means of a Zernike decomposition method

[S,C]=ZD[W _(L)(x,y)],  (19)

where S and C represent the sphere and cylinder power terms,respectively, and ZD represents a Zernike decomposition operator. Thevertex corrected sphere S′ and cylinder C′ power terms can be obtainedfrom the following formulae

$\begin{matrix}{{S^{\prime} = \frac{S}{1 + {dS}}},} & (20) \\{C^{\prime} = {\frac{S + C}{1 + {d\left( {S + C} \right)}} - {S^{\prime}.}}} & (21)\end{matrix}$

The vertex corrected wavefront can then be obtained by adding theuncorrected high order portion of the original wavefront with theZernike expansion operator applied to the corrected sphere S′ andcylinder C′ as

W′(x,y)=ZE(S′,C′)+W _(H)(x,y),  (22)

where ZE stands for a Zernike expansion operator.

b. Varying HOA

This algorithm segments the wavefront measurement into multipleportions, and is designed to have each portion of the correctedwavefront focused on or toward the focal point of the optical system,regardless of the wavefront shape. Thus, the local slope of each portionof the wavefront measurement can be scaled by a factor of f/(f+d) wheref represents the local focal length and d represents the vertexdistance. By applying the following algorithms, it is possible to obtainthe vertex corrected wavefront:

1. Calculate x- and y-gradient by the following algorithm:

-   -   Along the x axis:

∂W/∂x=[W(i,j+1)−W(i,j)]/dx if [i,j] lands on left edge  a.

∂W/∂x=[W(i,j)−W(i,j−1)]/dx if [i,j] lands on right edge  b.

∂W/∂x=[W(i,j+1)−W(i,j−1)]/2dx otherwise within pupil  c.

-   -   Along the y axis:

∂W/∂y=[W(i,j)−W(i+1,j)]/dy if [i,j] lands on upper edge  d.

∂W/∂y=[W(i−1,j)−W(i,j)]/dy if [i,j] lands on lower edge  e.

∂W/∂y=[W(i−1,j)−W(i+1,j)]/2dy otherwise within pupil  f.

-   -   If [i,j] is outside the pupil, the data is not considered.

2. Calculate local curvature P using this algorithm:

-   -   a. Calculate ∂²W/∂x², ∂²W/∂y² and ∂²W/∂x∂y from ∂W/∂x and ∂W/∂y        using algorithm 1.

$\frac{\partial^{2}W}{\partial r^{2}} = {{\frac{x^{2}}{x^{2} + y^{2}}\frac{\partial^{2}W}{\partial x^{2}}} + {\frac{2{xy}}{x^{2} + y^{2}}\frac{\partial^{2}W}{{\partial x}{\partial y}}} + {\frac{y^{2}}{x^{2} + y^{2}}\frac{\partial^{2}W}{\partial x^{2}}}}$

-   -   c. Calculate local curvature

$P = {\frac{1}{R^{2}}\frac{\partial^{2}W}{\partial r^{2}}}$

(R being pupil radius)

3. Scale the wavefront local curvature with this algorithm:

$\frac{\partial W^{\prime}}{\partial x} = {\frac{1}{1 + {Pd}}\frac{\partial W}{\partial x}}$$\frac{\partial W^{\prime}}{\partial y} = {\frac{1}{1 + {Pd}}\frac{\partial W}{\partial y}}$

4. Reconstruct the corrected wavefront W′(x,y) with this algorithm:

-   -   a. Calculate Fourier transform of ∂W′/θx and ∂W′/∂y to get c_(u)        and c_(v), respectively.    -   b. Multiply u with c_(u) and v with c_(v) and divide by u²+v².    -   c. Inverse Fourier transform to get W′(x,y).    -   d. Calculate ∂W′/∂x and ∂W′/∂y using algorithm 1, adjusted with        the edge being the entire frame as oppose to pupil edge.    -   e. Replace ∂W′/∂x and ∂W′/∂y with values from step 3 within the        pupil, leave values outside pupil untouched.    -   f. Determine if a predefined criteria is met, or if a        predetermined number of iterations have been completed. If not,        go to step (a) and repeat through step (f).    -   g. Provide an estimate of W′(x,y).

A predefined criteria of step (f) could be, for example, the RMS errorof the reconstructed wavefront based on a comparison between W′_(i) andW′_(i-1), in the ith and (i−1)th iterations, respectively.Alternatively, other optical quality gauges may be used. In oneembodiment, the predetermined number of iterations in step (f) is 10. Asillustrated in the above algorithm, it is possible to derive a treatmentsurface wavefront map based on an iterative Fourier reconstructionalgorithm. Thus the entire algorithm, steps 1 to 4, uses both Fourierreconstruction (step 4) and local slope scaling (step 3).

The theory behind Fourier reconstruction can be described as follows.Suppose wavefront W(x,y) is expanded into Fourier series as

W(x,y)=∫∫c(u,v)exp[i2π(ux+vy)]dudv,  (23)

where c(u,v) is the expansion coefficient. Taking partial derivative tox and y, respectively in the above equation, provides

$\begin{matrix}\left\{ \begin{matrix}{\frac{\partial{W\left( {x,y} \right)}}{\partial x} = {{2\pi}{\int{\int{{{uc}\left( {u,v} \right)}{\exp \left\lbrack {{2\pi}\left( {{ux} + {vy}} \right)} \right\rbrack}{u}{v}}}}}} \\{\frac{\partial{W\left( {x,y} \right)}}{\partial y} = {{2\pi}{\int{\int{{{vc}\left( {u,v} \right)}{\exp \left\lbrack {{2\pi}\left( {{ux} + {vy}} \right)} \right\rbrack}{u}{v}}}}}}\end{matrix} \right. & (24)\end{matrix}$

Denoting c_(u) to be the Fourier transform of x-derivative of W(x,y) andc_(v) to be the Fourier transform of y-derivative of W(x,y), provides

$\begin{matrix}\left\{ \begin{matrix}{\frac{\partial{W\left( {x,y} \right)}}{\partial x} = {\int{\int{{c_{u}\left( {u,v} \right)}{\exp \left\lbrack {\; 2\; {\pi \left( {{ux} + {vy}} \right)}} \right\rbrack}{u}{v}}}}} \\{\frac{\partial{W\left( {x,y} \right)}}{\partial y} = {\int{\int{{c_{v}\left( {u,v} \right)}{\exp\left\lbrack {\; 2\; {\pi \left( {{ux} + {vy}} \right)}{u}{v}} \right.}}}}}\end{matrix} \right. & (25)\end{matrix}$

Comparing these two sets of equations, provides

$\begin{matrix}\left\{ \begin{matrix}{{c_{u}\left( {u,v} \right)} = {\; 2\; \pi \; {{uc}\left( {u,v} \right)}}} \\{{c_{v}\left( {u,v} \right)} = {\; 2\; \pi \; {{vc}\left( {u,v} \right)}}}\end{matrix} \right. & (26)\end{matrix}$

Combining these two equations with u multiplied in both sides of thefirst equation and v multiplied in both sides of the second equation,provides

uc _(u)(u,v)+vc _(v)(u,v)=i2π(u ² v ²)c(u,v).  (27)

Therefore, the Fourier transform of wavefront can be obtained as

$\begin{matrix}\begin{matrix}{{c\left( {u,v} \right)} = {- \frac{\left\lbrack {{{uc}_{u}\left( {u,v} \right)} + {{vc}_{v}\left( {u,v} \right)}} \right\rbrack}{2\; {\pi \left( {u^{2} + v^{2}} \right)}}}} \\{= {- {\frac{}{2\; {\pi \left( {u^{2} + v^{2}} \right)}}\begin{bmatrix}{u{\int{\int{\frac{\partial{W\left( {x,y} \right)}}{\partial x}{\exp\left\lbrack {{{- }\; 2\; {\pi\left( {{ux} + {vy}} \right\rbrack}} +} \right.}}}}} \\{v{\int{\int{\frac{\partial{W\left( {x,y} \right)}}{\partial y}{\exp \left\lbrack {{- }\; 2\; {\pi \left( {{ux} + {vy}} \right)}} \right\rbrack}}}}}\end{bmatrix}}}}\end{matrix} & (28)\end{matrix}$

Hence, taking an inverse Fourier transform, it is possible to obtain thewavefront as

W(x,y)=∫∫c(u,v)exp[i2π(ux+vy)]dudv.  (29)

III. Refractive Treatment Shape

Once a treatment surface aberration has been derived by a method asdescribed above, it is possible to determine a prescription or arefractive treatment shape based on the treatment surface aberration.For example, a prescription can be derived for ameliorating a visioncondition in an eye of a patient. A refractive treatment shape can bedetermined based on the treatment surface aberration of the eye, and arefractive treatment shape can be embodied in any of a variety ofcorrective optical devices or procedures, including refractive lasersurgery, spectacles, contact lenses, intraocular lenses, and the like.

IV. Example Evaluating Classical Formulas and New Algorithms

In some embodiments, it is useful to evaluate the convergence of Fourierreconstruction used in the vertex correction algorithms discussed above.Such approaches are discussed in commonly owned patent application Ser.No. 10/601,048 filed Jun. 20, 2003, the entirety of which is herebyincorporated by reference. It is also useful to evaluate the accuracy ofthe varying high order aberration algorithm as compared to the classicalformulas discussed above (i.e. sphere, sphere and cylinder). Forexample, one test is to show the comparison between the algorithmicapproaches and the traditional approaches for myopic, hyperopic, andastigmatism cases. FIG. 7 shows the comparison of vertex correctedsphere and cylinder using the varying high order aberration algorithmdescribed above as compared to classical formulas (i.e. sphere, sphereand cylinder) for (a) hyperopia +3D; (b) myopia −3D; (c) astigmatism−2DS/−1.5DC. It is quite clear that the results are very good. Goodresults can be shown by a small error. For example, if the difference isless than 0.05 D, or smaller than 2.5%, it can generally be consideredgood. For pure sphere cases (e.g. myopia and hyperopia), the error canbe larger, due to coarse sampling of wavefront data in the calculation.

For high order aberrations, it has been shown with two examples (i.e.coma, spherical aberrations) in theory that the vertex correctedwavefront follows the power relationship given by the classical formulaof vertex correction. FIG. 8 shows wavefront surface plots of apre-vertex correction (left panel) and post-vertex correction (rightpanel) corresponding to a 12.5 mm vertex correction as accomplished by avarying high order aberration algorithm.

In terms of the efficiency of a varying high order aberration algorithm,the following table shows the running time taken for such a vertexcorrection algorithm with respect to the number of iterations taken inthe Fourier reconstruction, corresponding to step 4 of the algorithm, ina 1.13 GHz laptop computer. With 10 iterations, the algorithm can takemore than 2 seconds in real time, as shown in Table 1. Fortunately, thisvertex correction may only be needed when a treatment table isgenerated, which in itself may take minutes. Treatment tables are filesthat can store commands for a laser to deliver individual laser pulses,in the context of a laser ablation treatment. For example, the commandscan be for laser pulse duration and size.

TABLE 1 Iterations 1 2 5 10 20 50 200 Time (s) 0.340 0.521 1.231 2.3034.256 10.40 41.34

Thus in one embodiment, as part of the algorithm, Fourier reconstructioncan require about 10 iterations to achieve planned results given by100-micron sampling rate.

V. Wavefront Propagation from One Plane to Another

In some embodiments, the present invention provides treatment techniquesfor addressing high order aberrations, including algorithms that can beused to create treatment tables for custom ablation profiles. Thesetechniques can involve the determination of the expected target ablationprofile on the curved corneal surface, when the wavefront map is knownon the exit pupil plane, which is typically about 3.5 mm below thevertex of the corneal surface. In some cases, embodiments consider thewavefront to propagate as a whole, and do not involve addressing piecesof the wavefront separately. Algorithm embodiments can be validated witha classical vertex correction formula and a wavefront measuringexperiment. These high order aberration treatment techniques can beuseful in refining or determining shapes such as those used inophthalmic laser surgery. These techniques can be based on an opticalsystem, a human eye, with total ocular aberrations represented by awavefront on an exit pupil plane. In ray tracing terminology, this caninvolve a converging beam with potential aberration deviating from aperfect eye. These deviations can be modeled as wavefront aberrationsand be treated using the geometrical theory of aberrations.

A. Formulation of the Wavefront Propagation

As a wavefront propagates, it can be considered as many rays propagatingat different directions as determined by the norm of the local wavefrontsurface. According to the Huygens principle, the new wavefront is theenvelope of the spherical wavelets emanating from each point of theoriginal wavefront. FIG. 8A shows examples of a diverging wavefront anda wavefront that includes a spherical aberration before and afterpropagation according to embodiments of the present invention. Note thatthe wavefront boundary after propagation has been conformed to aslightly smaller area at the edge as the diffraction effect at the edgeis not a concern. In this section, a mathematical formulation is givenfor such a treatment. Examples of wavefront propagation according to theHuygens principle for FIG. 8A(a) a diverging defocus, and FIG. 8A(b) aspherical aberration are described.

1. Calculation of the Direction Factor

FIG. 8B shows the geometry of a myopic correction as the originalwavefront W (r, θ) with a radius R propagates a distance d from right toleft to become the new wavefront W′(r′, θ′) with a new radius R′. Thus,a geometry for a myopic wavefront with a radius R propagated a distanced from a less anterior plane towards a more anterior plane to a newwavefront with a radius R′ is provided. The reference plane for theoriginal wavefront is S and that for the new wavefront is S′. Thedirection at point T is determined by the angle between the norm of thewavefront at point T and the norm of the reference plane S, or the angleψ. This angle can be calculated from the radial slope of the wavefrontas

$\begin{matrix}{{\cos \; \psi} = {\frac{1}{\sqrt{1 + \left\lbrack \frac{\partial{W\left( {x,y} \right)}}{\partial x} \right\rbrack^{2} + \left\lbrack \frac{\partial{W\left( {x,y} \right)}}{\partial y} \right\rbrack^{2}}}.}} & ({A6})\end{matrix}$

For most applications, the wavefront slope is much smaller than 1. Forexample, even for a −10 D eye with a 6 mm pupil size, the maximum slopeis only 0.03, and its square is 0.0009. Hence, Eq. (A6) can beapproximated with a binomial expansion as

$\begin{matrix}{{\cos \; \psi} = {1 - {\frac{1}{2}\left\lbrack \frac{\partial{W\left( {x,y} \right)}}{\partial x} \right\rbrack}^{2} - {1{{\frac{1}{2}\left\lbrack \frac{\partial{W\left( {x,y} \right)}}{\partial y} \right\rbrack}^{2}.}}}} & ({A7})\end{matrix}$

Because Zernike polynomials use variables within a unit circle, the newvariables (ρ, θ) in polar coordinates and (u, v) in Cartesiancoordinates can be introduced in such a way that ρ=r/R and u=x/R, v=y/Rso that Eq. (A7) can be written as

$\begin{matrix}{{{\cos \; \psi} = {1 - {\frac{1}{2R^{2}}{a\left( {u,v} \right)}}}},} & ({A8})\end{matrix}$

where the direction factor a(u, v) can be written as

$\begin{matrix}{{a\left( {u,v} \right)} = {\left\lbrack \frac{\partial{W\left( {u,v} \right)}}{\partial u} \right\rbrack^{2} + {\left\lbrack \frac{\partial{W\left( {u,v} \right)}}{\partial v} \right\rbrack^{2}.}}} & ({A9})\end{matrix}$

From FIG. 8B, we have W=Q′Q and W′=T′T. In addition, d=TQ=T′P′. Letd′=T′Q′, we obtain

d′=d cos ψ.  (A10)

In addition, we have T′Q=T′T+TQ=W′+d=T′Q′+Q′Q=d′+W. Therefore,

$\begin{matrix}{{W^{\prime} - W} = {{d^{\prime} - d} = {{d\left( {{\cos \; \psi} - 1} \right)} = {{- \frac{d}{2R^{2}}}{{a\left( {u,v} \right)}.}}}}} & ({A11})\end{matrix}$

The magnitude of the difference in wavefronts can be represented asW′-W. The propagation distance can be represented as d. The directionfactor can be represented as a(u,v). R can represent the radius of theoriginal wavefront. Thus, the magnitude of the difference in wavefrontscan be proportional to the propagation distance. The magnitude of thedifference in wavefronts can also be proportional to the directionfactor. Similarly, the magnitude of the difference in wavefronts can beinversely proportional to the square of the wavefront radius or someother boundary dimension (for example the semi-major or semi-minor axisof an ellipse). FIG. 8C shows a geometry for a hyperopic wavefront witha radius R propagated a distance d from a less anterior plane towards amore anterior plane to a new wavefront with a radius R′. Because W′>0,W>0 and W′<W for the myopic case, Eq. (A11) is appropriate forrepresenting the propagation of a myopic wavefront. Similarly, for thehyperopic case, as shown in FIG. 8C, we have W=Q′Q and W′=T′T;d=T′Q′=P′Q; and d′=TQ. Therefore, we get T′Q=T′T+TQ=W′+d′=T′Q′+Q′Q=d+W.However, since W′<0, W<0 and |W′|>|W| for a hyperopic case, a negativesign needs to be used for W′ and W, or W+d′=d−W. This gives us

$\begin{matrix}{{W^{\prime} - W} = {{d^{\prime} - d} = {{d\left( {{\cos \; \psi} - 1} \right)} = {{- \frac{d}{2\; R^{2}}}{{a\left( {u,v} \right)}.}}}}} & ({A12})\end{matrix}$

which is identical to Eq. (A11). Therefore, Eq. (A11) can be used torepresent the propagation of any wavefront. Note that although themagnitude of the propagated wavefront is given by Eq. (A11), it isexpressed in the new coordinates (ρ′, θ′), not the original coordinates(ρ, θ), as discussed elsewhere herein. Again, following the previousconvention, d is negative when a wavefront is propagated from a moreanterior plane to a less anterior plane.

2. Calculation of the Propagated Zernike Coefficients and WavefrontBoundary

As further discussed elsewhere herein, in addition to calculation of thedirection factor, the formulation of the wavefront propagation caninclude calculation of the propagated Zernike coefficients and wavefrontboundary. Calculation of the propagated Zernike coefficients can involvethe use of Taylor monomials. Calculation of the propagated wavefrontboundary can involve a boundary factor.

3. Wavefront Propagation of a Converging Beam

In some embodiments, it can be helpful to consider a known wavefront onthe exit pupil plane for a myopic eye, represented by W(r,θ) in polarcoordinates where the optical path length with respect to the referencesphere S is given by Q′Q, as shown in FIG. 9. The propagation of the raythrough point Q is normal to the wavefront surface W at point Q. Whenthis ray travels a distance of d, the wavefront becomes W′(r′,θ′),because both the magnitude and coordinate system change. The ray may notnecessarily travel on the xz plane, although it is shown in FIG. 9 assuch. No matter which direction it goes, the ray can still be normal tothe new wavefront W′ and the new optical path length can be representedby T′T, with respect to the new reference sphere S′. If the anglebetween the normal of the wavefront W at Q and the normal of thereference sphere S at Q′ is φ, we have

d′=d cos φ.  (30)

FIG. 9 depicts a coordinate system for one point in the originalwavefront W(r,θ) with respect to the reference sphere S propagating toanother point in the new wavefront W′(r′,θ′) with respect to the newreference sphere S′ for a myopic eye, according to some embodiments.From FIG. 9, we obtain

Q′T=W′+d′=W+d,  (31)

where d′ is the distance of Q′T′ and d is the distance of QT and ofQ′A′. From Eqs. (30) and (31), we obtain

W′=W−d(cos⁻¹φ−1).  (32)

Because the slope of the wavefront can be related to the angle φ by

$\begin{matrix}{{{\cos \; \varphi} = {\frac{1}{\sqrt{1 + \left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2}}} \approx \frac{1}{1 + {\frac{1}{2}\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2}}}}},} & (33)\end{matrix}$

and the wavefront slope can be much smaller than 1 in some applications,from Eqs. (32) and (33) we obtain

$\begin{matrix}{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{W\left( {r,\theta} \right)} - {\frac{d}{2}\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2}}}} & (34)\end{matrix}$

Equation (34) indicates that the propagated wavefront can become smallerfor a myopic eye.

FIG. 10 shows a coordinate system for one point in the originalwavefront W(r,θ) with respect to the reference sphere S propagating toanother point in the new wavefront W′(r′,θ′) with respect to the newreference sphere S′ for a hyperopic eye, according to some embodiments.For a hyperopic eye, as shown in FIG. 10, the following relation can beobtained

QT'=W′+d=W+d′,  (35)

where d′ is the same as in Eq. (30). With similar processing aspreviously described, we obtain

$\begin{matrix}{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{W\left( {r,\theta} \right)} + {\frac{d}{2}{\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2}.}}}} & (36)\end{matrix}$

Equation (36) indicates that for a hyperopic eye, the propagatedwavefront can become larger. To simplify the calculation, Eq. (34) canbe used for both myopic and hyperopic cases, where −d is used instead ofd for hyperopic case. Similarly, for myopic case, if the propagatedplane is before the exit pupil plane, such as at spectacle plane, then−d can be used instead of d. For hyperopic eyes where the propagatedplane is before the exit pupil plane, d can be used instead of −d. NoteEqs. (34) and (36) can be independent of the focal length of theoriginal converging beam. In human eye case, that means they can beindependent of the power of the eye.

The approximation in Eq. (33) can be very small. Assuming a standard eyewith total power of 60 D, a 6 mm pupil, the wavefront slope fordifferent amounts of refractive power and the error in Eq. (33) areprovided in Table 2.

TABLE 2 Refractive power (D) Wavefront slope Error 1 0.0029889.96403E−12 2 0.005976 1.59468E−10 3 0.008965 8.07560E−10 4 0.0119552.55313E−09 5 0.014945 6.23540E−09 6 0.017936 1.29345E−08 7 0.0209282.39719E−08 8 0.023920 4.09115E−08 9 0.026914 6.55600E−08 10 0.0299089.99676E−08

B. Wavefront Propagation of a Parallel Beam

In some embodiments, it is possible to consider a parallel beam insteadof a converging beam, as shown in FIG. 11, which provides a coordinatesystem for one point in the original wavefront W(r,θ) with respect tothe reference plane S propagating to another point in the new wavefrontW′(r′,θ′) with respect to the new reference plane S′ for a myopic eye.Using a similar approach, we have

QT′=d−W′=d′−W,  (37)

which results in the following formula

$\begin{matrix}{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{W\left( {r,\theta} \right)} - {\frac{d}{2}{\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2}.}}}} & (38)\end{matrix}$

FIG. 12 provides a coordinate system for one point in the originalwavefront W(r,θ) with respect to the reference plane S propagating toanother point in the new wavefront W′(r′,θ′) with respect to the newreference plane S′ for a hyperopic eye, according to some embodiments.Similarly, for hyperopic eye, as shown in FIG. 12, we can obtain

$\begin{matrix}{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{W\left( {r,\theta} \right)} + {\frac{d}{2}{\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2}.}}}} & (39)\end{matrix}$

According to some embodiments of the present invention, the formulae forparallel beams can be the same as those for converging beams. Asdiscussed elsewhere herein, in some cases use of a parallel beam cansimplify the treatment of wavefront propagation.

Even though the cornea has a strong curvature, when using a wavefrontdevice it is possible to assume parallel propagation of a wavefront.This may be true even where a patient has myopia or hyperopia, forexample.

C. Wavefront Boundary or Size

As noted previously, the wave-front boundary can change when thewavefront propagates. For example, for low order spherocylindricalerror, a circular wavefront becomes an elliptical wavefront when itpropagates, as shown in FIG. 12A. A low order wavefront (c₂ ⁻²=1 μm, c₂⁰=3 μm, c₂ ²=2 μm, and R=3 mm) propagated to become an ellipticalwavefront, is shown in FIG. 12A, where FIG. 12A(a) represents d=35 mm,FIG. 12A(b) represents d=125 mm, and FIG. 12A(c) represents d=120 mm.The aspect ratios of the ellipses are 0.9914, 0.9691, and 0.6638,respectively. In Atchison, D. A. et al., J. Opt. Soc. Am. A 20, 965-973(2003), a subject was investigated that deals with the boundary changefrom a circular pupil to an elliptical pupil when it is off-axis.Similarly, for a coma wavefront, it also becomes elliptical. For asecondary astigmatism, it becomes bi-elliptical, or a fourfold symmetry.The instant discussion describes an exemplary approach for calculatingthe propagated wavefront boundary.

From FIG. 8B, the relationship between the original wavefront radius Rand the propagated wavefront radius R′ can be determined from the simplegeometry as

$\begin{matrix}{{{R^{\prime} - R} = {{{- d}\; \tan \; \psi} = {{- d}\frac{b}{R^{\prime}}}}},} & ({A17})\end{matrix}$

where the boundary factor b can be written as

$\begin{matrix}\begin{matrix}{b = {\left\lbrack {\left( \frac{\partial{W\left( {u,v} \right)}}{\partial u} \right)^{2} + \left( \frac{\partial{W\left( {u,v} \right)}}{\partial v} \right)^{2}} \right\rbrack^{1/2}_{\rho = 1}}} \\{= {\sqrt{a\left( {u,v} \right)}_{\rho = 1}.}}\end{matrix} & ({A18})\end{matrix}$

Therefore, the boundary factor b is the square root of the directionfactor at the boundary of the original wavefront (i.e., ρ=√{square rootover (u²+v²)}=1). Similarly, from FIG. 8C, we find a similarrelationship as

$\begin{matrix}{{R^{\prime} - R} = {{d\; \tan \; \psi} = {{- d}{\frac{b}{R}.}}}} & ({A19})\end{matrix}$

For a converging wavefront, the wavefront radius becomes larger when itpropagates backwards, but the defocus coefficient is negative. For adiverging wavefront, the wavefront radius becomes smaller, but thedefocus coefficient is positive. Therefore, Eq. (A17) can be used forcases when the signs of Zernike coefficients are correctly applied.

In vision applications, the low order spherocylindrical error istypically much larger than the high order ocular aberrations. Therefore,the influence of the wavefront propagation on the new wavefront boundaryfor the low order aberrations may be much more significant than that forthe high order aberrations.

Dai, G.-m. J. Opt. Soc. Am. A 23:1657-1666 (2006) provides a method todefine the Zernike polynomials as

Z _(i)(ρ,θ)=

|₂ ^(m)|(ρ)Θ^(m)(θ),  (A20)

where n and m denote the radial degree and the azimuthal frequency,respectively, the radial polynomials are defined as

n  m   ( ρ ) = ∑ s = 0 ( n -  m  ) / 2  ( - 1 ) s  n + 1  ( n -s ) !  ρ n - 2  s s !  [ ( n + m ) / 2 - s ] !  [ ( n - m ) / 2 - s] ! , ( A21 )

and the triangular functions as

$\begin{matrix}{{\Theta^{m}(\theta)} = \left\{ \begin{matrix}{\sqrt{2}\cos {m}\theta} & \left( {m > 0} \right) \\1 & \left( {m = 0} \right) \\{\sqrt{2}\sin {m}\theta} & {\left( {m < 0} \right).}\end{matrix} \right.} & ({A22})\end{matrix}$

If we represent the direction factor a(u, v) with Zernike polynomials,we may write

$\begin{matrix}{{a\left( {u,v} \right)} = {{a\left( {p,\theta} \right)} = {\sum\limits_{i = 1}^{J}{g_{i}{{Z_{i}\left( {\rho,\theta} \right)}.}}}}} & ({A23})\end{matrix}$

We can separate the radially symmetric terms and the radially asymmetricpairs of terms as

a  ( ρ , θ ) =  ∑ n , m  n  m   ( ρ )  2  ( g n -  m   sin  m   θ + g n  m   cos   m   θ ) +  ∑ n  g n 0  n 0  ( ρ )=  ∑ n , m m ≠ 0 m ≠ 0   n  m   ( ρ )  2  ( g n -  m  ) 2 + (g n  m  ) 2  cos   m   θ - φ ) +  ∑ n  g n 0  n 0  ( ρ ) , (A24 )

where the angle of |m|-symmetry φ can be expressed as

$\begin{matrix}{\varphi = {\frac{1}{m}{{\tan^{- 1}\left( \frac{c_{n}^{- {m}}}{c_{n}^{m}} \right)}.}}} & ({A25})\end{matrix}$

It can be shown (Born, M. et al., Principles of Optics, 7th ed.(Cambridge University Press, 1999)) that

_(n) ^(m)(1)=√{square root over (n−1)}.  (A26)

Note that as described herein, in one embodiment an exemplary definitionof Zernike radial polynomials can differ from the definition in Born etal. by a factor of √{square root over (n+1)}.

With the use of Eq. (A24), the boundary factor b can be written as

$\begin{matrix}{b^{2} = {{\sum\limits_{n}^{\;}{\sqrt{n + 1}g_{n}^{0}}} + {\sum\limits_{\underset{m \neq 0}{n,m}}^{\;}{\sqrt{2\left( {n + 1} \right)}\sqrt{\left( g_{n}^{- {m}} \right)^{2} + \left( g_{n}^{m} \right)^{2}}\cos {m}{\left( {\theta - \varphi} \right).}}}}} & ({A27})\end{matrix}$

The average of the triangular function of Eq. (A27) is zero. Therefore,as first order approximation, Eq. (A27) can be written as

$\begin{matrix}{b^{2} = {\sum\limits_{n}^{\;}{\sqrt{n + 1}{g_{n}^{0}.}}}} & ({A28})\end{matrix}$

Hence, the final propagated wavefront has a radius R′ as

$\begin{matrix}{R^{\prime} = {{R\left( {1 - {d\frac{b}{R^{2}}}} \right)}.}} & ({A29})\end{matrix}$

As described previously, it is possible to use the reference spheres orplanes S and S′ to determine the optical path length and to assume thewavefront comes from either a converging beam or a parallel beam. Insome cases where the wavefront comes from a converging beam, the neteffect of the beam size due to the aberrations can be complicated. Ifthe wavefront comes from a parallel beam, however, the net effect of thebeam size can be simpler because the beam is parallel and so does notchange the beam size. The effective beam size can thus solely bedetermined by the aberration. In some embodiments, with spectaclecorrection, a myope may see things smaller and a hyperope may see thingslarger after the correction. This can be described as minification andmagnification effects. From FIG. 11, the right triangle QA′T indicates

$\begin{matrix}{{r - r^{\prime}} = {{{d\; \sin \; \varphi} \approx {d\; \tan \; \varphi}} = {d{\frac{\partial{W\left( {r,\theta} \right)}}{\partial r}.}}}} & (40)\end{matrix}$

Note the approximation in Eq. (40) is good, as the error (tan φ−sin φ)is small, as can be seen from Table 3.

TABLE 3 Refractive power (D) Wavefront slope Error 1 0.0029881.33386E−08 2 0.005976 1.06706E−07 3 0.008965 3.60242E−07 4 0.0119558.54225E−07 5 0.014945 1.66873E−06 6 0.017936 2.88431E−06 7 0.0209284.58153E−06 8 0.023920 6.84017E−06 9 0.026914 9.74247E−06 10 0.0299081.33672E−05

Equation (40) indicates that for a myopic eye, the wavefront size canbecome smaller in the propagation direction and the change of the sizecan be proportional to the propagation distance and the wavefront slope.Similarly, for a hyperopic eye, we have

$\begin{matrix}{{r - r^{\prime}} = {{{d\; \sin \; \varphi} \approx {d\; \tan \; \varphi}} = {d{\frac{\partial{W\left( {r,\theta} \right)}}{\partial r}.}}}} & (41)\end{matrix}$

Now, consider a balanced defocus term, represented with Zernikepolynomials as

$\begin{matrix}{{{W\left( {r,\theta} \right)} = {{c_{2}^{0}{Z_{2}^{0}\left( {r,\theta} \right)}} = {\sqrt{3}{c_{2}^{0}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack}}}},} & (42)\end{matrix}$

where R is the radius of the aperture on the exit pupil. From Eq. (42)we get

$\begin{matrix}{{\frac{\partial{W\left( {r,\theta} \right)}}{\partial r} = {4\sqrt{3}c_{2}^{0}\frac{r}{R^{2}}}},} & (43)\end{matrix}$

so at the periphery of the wavefront, i.e., r=R, we have

$\begin{matrix}{{\frac{\partial{W\left( {r,\theta} \right)}}{\partial r}_{r = R}} = {4\sqrt{3}c_{2}^{0}{\frac{1}{R}.}}} & (44)\end{matrix}$

Substituting Eq. (44) into (40) for the periphery of the wavefront, wehave

$\begin{matrix}{R^{\prime} = {{R - {d\frac{\partial{W\left( {r,\theta} \right)}}{\partial r}}} = {{R\left( {1 - {d\frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}} \right)}.}}} & (45)\end{matrix}$

Similarly, for a hyperopic eye, we have

$\begin{matrix}{R^{\prime} = {{R + {d\frac{\partial{W\left( {r,\theta} \right)}}{\partial r}}} = {{R\left( {1 + {d\frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}} \right)}.}}} & (46)\end{matrix}$

Again, we can just use Eq. (45) for both cases where for a hyperopiceye, −d is used instead of d.

D. Sphere

It is possible to show that the combination of Eqs. (38) and (41) for apure defocus can give the same formula as the classical vertexcorrection formula

$\begin{matrix}{{S^{\prime} = \frac{S}{1 + {Sd}}},} & (47)\end{matrix}$

where S′ and S are in diopters and d in meters. Substituting Eq. (42)into (38), we have

$\begin{matrix}{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{c_{2}^{0}\left\lbrack {1 - {d\frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}} \right\rbrack}{{\sqrt{3}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack}.}}} & (48)\end{matrix}$

From Eqs. (37), (39), and (41), we know

$\begin{matrix}{\frac{r^{\prime}}{R^{\prime}} = {\frac{r}{R}.}} & (49)\end{matrix}$

Substituting Eq. (49) into (48), we obtain

$\begin{matrix}{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{{c_{2}^{0}\left\lbrack {1 - {d\frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}} \right\rbrack}{\sqrt{3}\left\lbrack {{2\left( \frac{r^{\prime}}{R^{\prime}} \right)^{2}} - 1} \right\rbrack}} - {\frac{12{d\left( c_{2}^{0} \right)}^{2}}{R^{2}}.}}} & (50)\end{matrix}$

If we use the normalized radial variable

$\begin{matrix}{{\rho = \frac{r}{R}},{\rho^{\prime} = {\frac{r^{\prime}}{R^{\prime}}.}}} & (51)\end{matrix}$

we can rewrite Eqs. (38) and (50) as

$\begin{matrix}{{{W\left( {\rho,\theta} \right)} = {{c_{2}^{0}{Z_{2}^{0}\left( {\rho,\theta} \right)}} = {\sqrt{3}{c_{2}^{0}\left( {{2\rho^{2}} - 1} \right)}}}},} & (52) \\{{{W^{\prime}\left( {\rho^{\prime},\theta^{\prime}} \right)} = {{\sqrt{3}{b_{2}^{0}\left( {{2\rho^{\prime 2}} - 1} \right)}} + b_{0}^{0}}},{where}} & (53) \\{{b_{2}^{0} = {c_{2}^{0}\left\lbrack {1 - {d\frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}} \right\rbrack}},} & (54) \\{b_{0}^{0} = {- {\frac{12{d\left( c_{2}^{0} \right)}^{2}}{R^{2}}.}}} & (55)\end{matrix}$

In some cases, the induced piston term does not have any significance onimage quality and can be ignored. In addition, it may be very small.From the definition of wavefront refractions

$\begin{matrix}{{S = {- \frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}},} & (56) \\{S^{\prime} = {- {\frac{4\sqrt{3}b_{2}^{0}}{R^{\prime 2}}.}}} & (57)\end{matrix}$

Substituting Eq. (55) into Eq. (56) and performing certain arithmeticoperations, we get

$\begin{matrix}{S^{\prime} = {\frac{S}{1 + {Sd}}.}} & (58)\end{matrix}$

E. Sphere and Cylinder

It is possible to show that the combination of Eqs. (34) and (41) for apure defocus can give the same formula as the classical vertexcorrection formula Eq. (58) and

$\begin{matrix}{{{S^{\prime} + C^{\prime}} = \frac{S + C}{1 + {d\left( {S + C} \right)}}},} & (59)\end{matrix}$

where S′, C′, S, and S′ are in diopters and d in meters. Now thewavefront containing the sphere and cylinder can be represented byZernike polynomials Z₃, Z₄, and Z₅ as

$\begin{matrix}{{W\left( {r,\theta} \right)} = {{c_{2}^{- 2}\sqrt{6}\left( \frac{r}{R} \right)^{2}\sin \; 2\; \theta} + {c_{2}^{0}{\sqrt{3}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack}} + {c_{2}^{2}\sqrt{6}\left( \frac{r}{R} \right)^{2}\cos \; 2{\theta.}}}} & (60)\end{matrix}$

Taking a derivative of W(r,θ) with respect to x and to y, we obtain

$\begin{matrix}{{\frac{\partial{W\left( {r,\theta} \right)}}{\partial x} = {{2\sqrt{6}c_{2}^{- 2}\frac{1}{R}\left( \frac{r}{R} \right)\sin \; \theta} + {4\sqrt{3}c_{2}^{0}\frac{1}{R}\left( \frac{r}{R} \right)\cos \; \theta} + {2\sqrt{6}c_{2}^{2}\frac{1}{R}\left( \frac{r}{R} \right)\cos \; \theta}}},} & (61) \\{\frac{\partial{W\left( {r,\theta} \right)}}{\partial y} = {{2\sqrt{6}c_{2}^{- 2}\frac{1}{R}\left( \frac{r}{R} \right)\cos \; \theta} + {4\sqrt{3}c_{2}^{0}\frac{1}{R}\left( \frac{r}{R} \right)\sin \; \theta} - {2\sqrt{6}c_{2}^{2}\frac{1}{R}\left( \frac{r}{R} \right)\sin \; {\theta.}}}} & (62)\end{matrix}$

Therefore, we get Eq. (63), as follows.

$\begin{matrix}\begin{matrix}{\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial r} \right)^{2} = {\left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial x} \right)^{2} + \left( \frac{\partial{W\left( {r,\theta} \right)}}{\partial y} \right)^{2}}} \\{= {{\frac{16\sqrt{3}}{R^{2}}c_{2}^{- 2}c_{2}^{0}\sqrt{6}\left( \frac{r}{R} \right)^{2}\sin \; 2\; \theta} + \frac{4\sqrt{3}}{R^{2}}}} \\{{{\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack {\sqrt{3}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack}} +}} \\{{{\frac{16\sqrt{3}}{R^{2}}c_{2}^{2}c_{2}^{0}\sqrt{6}\left( \frac{r}{R} \right)^{2}\cos \; 2\; \theta} + \frac{12}{R^{2}}}} \\{{\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack.}} \\{= {\frac{4\sqrt{3}}{R^{2}}\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}} \\{{{\sqrt{3}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack} + \frac{12}{R^{2}}}} \\{{\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack +}} \\{{\frac{16\sqrt{3}}{R^{2}}\sqrt{\left( c_{2}^{- 2} \right)^{2} - \left( c_{2}^{2} \right)^{2}}}} \\{{c_{2}^{0}\sqrt{6}\left( \frac{r}{R} \right)^{2}\cos \; 2{\left( {\theta - \varphi} \right).}}}\end{matrix} & (63)\end{matrix}$

From Eq. (63), it may be shown that due to the propagation of asymmetricwavefront, the wavefront boundary can become elliptical, for thewavefront slopes at different meridian can be different. The maximum andminimum wavefront slopes can be obtained by setting cos 2(θ−φ) to 1 and−1 and r=R, or

$\begin{matrix}{{\frac{\partial{W\left( {r,\theta} \right)}}{\partial r}_{r = R}} = {{\frac{1}{R}\left\lbrack {{4\sqrt{3}c_{2}^{0}} \pm {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} - \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}.}} & (64)\end{matrix}$

So the maximum and minimum wavefront radius can be calculated from Eq.(45) as

$\begin{matrix}{R^{\prime} = {{R\left\lbrack {1 - {\frac{d}{R}\left( {{4\sqrt{3}c_{2}^{0}} \mp {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right)}} \right\rbrack}.}} & (65)\end{matrix}$

From Eqs. (37), (39), and (41), we know

$\begin{matrix}{\frac{r^{\prime}}{R^{\prime}} = {\frac{r}{R}.}} & (66)\end{matrix}$

With no rotation, we have

θ′=θ.  (67)

Substituting Eqs. (63)-(67) and (60) into (38), we have

$\begin{matrix}{{{W^{\prime}\left( {r^{\prime},\theta^{\prime}} \right)} = {{b_{2}^{- 2}\sqrt{6}\left( \frac{r^{\prime}}{R^{\prime}} \right)^{2}\sin \; 2\; \theta^{\prime}} + {b_{2}^{0}{\sqrt{3}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack}} + {b_{2}^{2}\sqrt{6}\left( \frac{r^{\prime}}{R^{\prime}} \right)^{2}\cos \; 2\; \theta^{\prime}} + b_{0}^{0}}},\mspace{20mu} {where}} & (68) \\{\mspace{79mu} {{b_{2}^{- 2} = {c_{2}^{- 2}\left\lbrack {1 - {d\frac{8\sqrt{3}c_{2}^{0}}{R^{2}}}} \right\rbrack}},}} & (69) \\{\mspace{79mu} {{b_{2}^{0} = {c_{2}^{0}\left\lbrack {1 - {d\frac{2\sqrt{3}\left\{ {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\}}{c_{2}^{0}R^{2}}}} \right\rbrack}},}} & (70) \\{\mspace{79mu} {{b_{2}^{2} = {c_{2}^{2}\left\lbrack {1 - {d\frac{8\sqrt{3}c_{2}^{0}}{R^{2}}}} \right\rbrack}},}} & (71) \\{\mspace{79mu} {b_{0}^{0} = {{\frac{12}{R^{2}}\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}.}}} & (72)\end{matrix}$

When using Eq. (49) for the wavefront as represented in Eq. (60), it canbe found that the new pupil radius may no longer be circular, because ofthe astigmatism term. For wavefront representation within a circulararea, the defocus term can be used, which results in the same formula asEq. (49) for the new pupil radius. In addition, the piston induced dueto the wavefront propagation may not have any imaging consequence. Usinga plus cylinder notation, we can rewrite Eqs. (60) and (68) as

$\begin{matrix}{{{W\left( {r,\theta} \right)} = {{c_{2}^{0}{\sqrt{3}\left\lbrack {{2\left( \frac{r}{R} \right)^{2}} - 1} \right\rbrack}} + {\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}\sqrt{6}\left( \frac{r}{R} \right)^{2}\cos \; 2\left( {\theta - \varphi} \right)}}},} & (73) \\{{{W\left( {r^{\prime},\theta^{\prime}} \right)} = {{b_{2}^{0}{\sqrt{3}\left\lbrack {{2\left( \frac{r^{\prime}}{R^{\prime}} \right)^{2}} - 1} \right\rbrack}} + {\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}\sqrt{6}\left( \frac{r^{\prime}}{R^{\prime}} \right)^{2}\cos \; 2\left( {\theta^{\prime} - \varphi^{\prime}} \right)}}},} & (74)\end{matrix}$

where φ and φ′ are the cylinder axis for W and W′, respectively, and aregiven by

$\begin{matrix}{{\varphi = {\frac{1}{2}{\tan^{- 1}\left( \frac{c_{2}^{- 2}}{c_{2}^{2}} \right)}}},} & (75) \\{\varphi^{\prime} = {{\frac{1}{2}{\tan^{- 1}\left( \frac{b_{2}^{- 2}}{b_{2}^{2}} \right)}} = {{\frac{1}{2}{\tan^{- 1}\left( \frac{c_{2}^{- 2}}{c_{2}^{2}} \right)}} = {\varphi.}}}} & (76)\end{matrix}$

Because the refraction is related to the wavefront curvature, the powercan be written as

$\begin{matrix}{{{P\left( {r,\theta} \right)} = {{- \frac{4\sqrt{3}c_{2}^{0}}{R^{2}}} - \frac{2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}{R^{2}} + {\frac{2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}{R^{2}}\left\lbrack {1 - {\cos \; 2\left( {\theta - \varphi} \right)}} \right\rbrack}}},} & (77) \\{{P\left( {r^{\prime},\theta^{\prime}} \right)} = {{- \frac{4\sqrt{3}b_{2}^{0}}{R^{\prime 2}}} - \frac{2\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}{R^{\prime 2}} + {{\frac{2\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}{R^{\prime 2}}\left\lbrack {1 - {\cos \; 2\left( {\theta^{\prime} - \varphi^{\prime}} \right)}} \right\rbrack}.}}} & (78)\end{matrix}$

Hence, the sphere and cylinder are

$\begin{matrix}{{C = \frac{4\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}{R^{2}}},} & (79) \\{{S = {\frac{4\sqrt{3}c_{2}^{0}}{R^{2}} - \frac{2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}{R^{2}}}},} & (80) \\{{C^{\prime} = \frac{4\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}{R^{\prime 2}}},} & (81) \\{S^{\prime} = {{- \frac{4\sqrt{3}b_{2}^{0}}{R^{\prime 2}}} - {\frac{2\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}{R^{\prime 2}}.}}} & (82)\end{matrix}$

When we consider the minimum power (or sphere only), the new wavefrontradius from Eq. (65) can be written as

$\begin{matrix}{R^{\prime} = {{R\left\lbrack {1 - {\frac{d}{R^{2}}\left( {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right)}} \right\rbrack} = {{R\left( {1 + {dS}} \right)}.}}} & (83)\end{matrix}$

Expanding Eq. (82) with some algebra, we get

$\begin{matrix}\begin{matrix}{S^{\prime} = {- {\frac{1}{R^{\prime 2}}\left\lbrack {{4\sqrt{3}b_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}}} \\{= {- {\frac{1}{R^{\prime 2}}\left\lbrack {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}} -} \right.}}} \\\left. {\frac{d}{R^{2}}\left( {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right)^{2}} \right\rbrack \\{= {- {\frac{1}{{R^{2}\left( {1 + {dS}} \right)}^{2}}\left\lbrack {{- R^{2}}{S\left( {1 + {dS}} \right)}} \right\rbrack}}} \\{= {\frac{S}{1 + {dS}}.}}\end{matrix} & (84)\end{matrix}$

Similarly, if we consider the maximum power (or sphere plus cylinder),we have the new wavefront radius from Eq. (65) as

$\begin{matrix}\begin{matrix}{R^{\prime} = {R\left\lbrack {1 - {\frac{d}{R^{2}}\left( {{4\sqrt{3}c_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right)}} \right\rbrack}} \\{= {{R\left( {1 + {d\left( {S + C} \right)}} \right\rbrack}.}}\end{matrix} & (85)\end{matrix}$

Expanding the sum of Eqs. (81) and (82), we obtain

$\begin{matrix}\begin{matrix}{{S^{\prime} + C^{\prime}} = {- {\frac{1}{R^{\prime 2}}\left\lbrack {{4\sqrt{3}b_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}}} \\{= {- {\frac{1}{R^{\prime 2}}\left\lbrack {{4\sqrt{3}c_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}} -} \right.}}} \\\left. {\frac{d}{R^{2}}\left( {{4\sqrt{3}c_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right)^{2}} \right\rbrack \\{= {{- \frac{1}{{R^{2}\left\lbrack {1 + {d\left( {S + C} \right)}} \right\rbrack}^{2}}}\left\{ {- {{R^{2}\left( {S + C} \right)}\left\lbrack {1 + {d\left( {S + C} \right)}} \right\rbrack}} \right\}}} \\{= {\frac{S + C}{1 + {d\left( {S + C} \right)}}.}}\end{matrix} & (86)\end{matrix}$

Hence Eqs. (58) and (59) can be proven.

F. General Analytical Representation

According to some implementation embodiments of the present invention,it may be helpful to consider the new wavefront radius that is half ofthe maximum and minimum. In this case, we consider the sphericalequivalent as the power, or the new defocus term, so the new wavefrontradius can be expressed as

$\begin{matrix}{R^{\prime} = {{R\left\lbrack {1 - {d\frac{4\sqrt{3}c_{2}^{0}}{R^{2}}}} \right\rbrack}.}} & (87)\end{matrix}$

For ocular aberrations, the defocus term is often the dominant term.Table 4 lists the error in percentage when Eq. (87) is used instead of astandard formula (Plus cylinder notation is used).

TABLE 4 S (D) C (D) d (mm) Percent Error −8 +6 3.5 ±1.07% 0 +6 3.5±1.04% +2 +4 3.5 ±0.69% −8 +6 16 ±5.22% +2 +4 16 ±3.01%

For a propagation distance of 3.5 mm, that is the distance from the exitpupil plane to the cornea plane, the error in terms of the pupil radiusat extreme cylinder case is barely one percent. To the spectacle plane,it becomes somewhat significant, which stands at 5% maximum. Therefore,for the purpose of wavefront propagation from exit pupil plane to thecorneal plane, Eq. (87) can be used as an approximation. For ocularaberrations, modal representation with the use of Zernike polynomials isoften used as

$\begin{matrix}{{{W\left( {r,\theta} \right)} = {{\sum\limits_{i = 1}^{J}{a_{i}{Z_{i}\left( {r,\theta} \right)}}} = {\sum\limits_{n = 1}^{N}{\sum\limits_{m = {- n}}^{n}{a_{n}^{m}{R_{n}^{m}(r)}{\Theta^{m}(\theta)}}}}}},} & (88)\end{matrix}$

where a_(i) is the coefficient of the ith polynomials, and the radialpolynomials

$\begin{matrix}{{R_{n}^{m}(r)} = {\sum\limits_{s = 0}^{{({n - m})}/2}{\frac{\left( {- 1} \right)^{s}\sqrt{n + 1}{\left( {n - s} \right)!}}{{{{{s!}\left\lbrack {{\left( {n + m} \right)/2} - s} \right\rbrack}!}\left\lbrack {{\left( {n - m} \right)/2} - s} \right\rbrack}!}\left( \frac{r}{R} \right)^{n - {2s}}}}} & (89)\end{matrix}$

and the triangle functions

$\begin{matrix}{{\Theta^{m}(\theta)} = \left\{ \begin{matrix}{\sqrt{2}\sin \; m\; \theta} & \left( {m < 0} \right) \\1 & \left( {m = 0} \right) \\{\sqrt{2}\cos \; m\; \theta} & \left( {m > 0} \right)\end{matrix} \right.} & (90)\end{matrix}$

In some cases, the complex representation of wavefront with Zernikepolynomials can make it difficult, if not impossible, to calculate thesquares of the wavefront slopes. On the other hand, the simple form ofTaylor monomials can make it easy to calculate the squares of thewavefront slopes for the calculation of the direction angle φ as

$\begin{matrix}{{\cos^{- 1}\varphi} = {1 + {\frac{1}{2R^{2}}\left( \frac{\partial{W\left( {x,y} \right)}}{\partial x} \right)^{2}} + {\frac{1}{2\; R^{2}}{\left( \frac{\partial{W\left( {x,y} \right)}}{\partial y} \right)^{2}.}}}} & (91)\end{matrix}$

If the wavefront is expanded into Taylor monomials, we have

$\begin{matrix}{{{W\left( {r,\theta} \right)} = {{\sum\limits_{i = 1}^{J}{\alpha_{i}{T_{i}\left( {r,\theta} \right)}}} = {\sum\limits_{p,q}{\frac{\alpha_{p}^{q}}{R^{p}}{T_{p}^{q}\left( {x,y} \right)}}}}},} & (92)\end{matrix}$

where Taylor monomials are defined as

T _(p) ^(q)(ρ,θ)=ρ^(p) cos^(q)θ sin^(p-q) θ=T _(p) ^(q)(x,y)=x ^(q) y^(p-q).  (93)

This simple form makes Eq. (91) representable by a linear combination ofTaylor monomials as Eq. (94) as follows.

$\begin{matrix}{{\cos^{- 1}\varphi} = {1 + {\frac{1}{2\; R^{2}}{\sum\limits_{p,q}{\sum\limits_{p^{\prime},q^{\prime}}{\alpha_{p}^{q}\alpha_{p^{\prime}}^{q^{\prime}}{qq}^{\prime}{T_{p + p^{\prime} - 2}^{q + q^{\prime} - 2}\left( {x,y} \right)}}}}} +}} \\{{\frac{1}{2R^{2}}{\sum\limits_{p,q}{\sum\limits_{p^{\prime},q^{\prime}}{\alpha_{p}^{q}{\alpha_{p^{\prime}}^{q^{\prime}}\left( {p - q} \right)}\left( {p^{\prime} - q^{\prime}} \right){T_{p + p^{\prime} - 2}^{q + q^{\prime}}\left( {x,y} \right)}}}}}} \\{= {1 + {\frac{1}{2\; R^{2}}{\sum\limits_{i = 1}^{J^{\prime}}{\beta_{i}{{T_{i}\left( {x,y} \right)}.}}}}}}\end{matrix}$

Embodiments of the present invention encompass various techniques thatinvolve the calculation of propagated Zernike coefficients. For example,it is possible to obtain an analytical expression of the directionfactor in terms of Zernike polynomials. It may be helpful to use Taylormonomials because a(u, v) can be obtained for a wavefront with Taylormonomials (Dai, G.-m. J. Opt. Soc. Am. A 23:2970-2971 (2006)) as

$\begin{matrix}{{\left\lbrack \frac{\partial{W\left( {u,v} \right)}}{\partial u} \right\rbrack^{2} = {\sum\limits_{p,q}{\sum\limits_{p^{\prime},q^{\prime}}{\alpha_{p}^{q}\alpha_{p^{\prime}}^{q^{\prime}}{qq}^{\prime}{T_{p + p^{\prime} - 2}^{q + q^{\prime} - 2}\left( {u,v} \right)}}}}},} & ({A13a}) \\{{\left\lbrack \frac{\partial{W\left( {u,v} \right)}}{\partial u} \right\rbrack^{2} = {\sum\limits_{p,q}{\sum\limits_{p^{\prime},q^{\prime}}{\alpha_{p}^{q}{\alpha_{p^{\prime}}^{q^{\prime}}\left( {p - q} \right)}\left( {p^{\prime} - q^{\prime}} \right){T_{p + p^{\prime} - 2}^{q + q^{\prime}}\left( {u,v} \right)}}}}},} & \left( {13b} \right)\end{matrix}$

where α_(p) ^(q) is the Taylor coefficient when the wavefront isexpanded into Talor monomials as

$\begin{matrix}\begin{matrix}{{W\left( {\rho,\theta} \right)} = {\sum\limits_{i = 1}^{J}{\alpha_{i}{T_{i}\left( {\rho,\theta} \right)}}}} \\{= {\sum\limits_{p,q}{\alpha_{p}^{q}\rho^{p}\cos^{q}\theta \; \sin^{p - q}\theta}}} \\{{= {\sum\limits_{p,q}{\alpha_{p}^{q}u^{q}v^{p - q}}}},}\end{matrix} & ({A14})\end{matrix}$

where J is the total number of Taylor monomials in the wavefrontexpansion. In Eq. (A14), we have used both the single index i and thedouble index (p, q) for Taylor monomials. Similarly to the indeces ofZernike polynomials, p is referred to as the radial order and q theazimuthal frequency (Dai, G.-m. J. Opt. Soc. Am. A 23:1657-1666 (2006).From Eqs. (A8), (A13) and (A14), we obtain

$\begin{matrix}\begin{matrix}{{a\left( {u,v} \right)} = {{\sum\limits_{p,q}{\sum\limits_{p^{\prime},q^{\prime}}{\alpha_{p}^{q}\alpha_{p^{\prime}}^{q^{\prime}}{qq}^{\prime}{T_{p + p^{\prime} - 2}^{q + q^{\prime} - 2}\left( {u,v} \right)}}}} +}} \\{{\sum\limits_{p,q}{\sum\limits_{p^{\prime},q^{\prime}}{\alpha_{p}^{q}{\alpha_{p^{\prime}}^{q^{\prime}}\left( {p - q} \right)}\left( {p^{\prime} - q^{\prime}} \right){T_{p + p^{\prime} - 2}^{q + q^{\prime} - 2}\left( {u,v} \right)}}}}} \\{{= {\sum\limits_{i = 1}^{J^{\prime}}{\beta_{i}{T_{i}\left( {u,v} \right)}}}},}\end{matrix} & ({A15})\end{matrix}$

where J′ is the new number of monomials to be affected by the wavefrontpropagation and β_(i) is the coefficient of the ith monomial after thepropagation. Table 5 shows the conversion table for α_(i) to β_(i) forJ=27 (6th order). It can be shown that the new radial degree n′=2n−2,where n is the original radial degree. For example, if J=27 (6th order),then J′=65 (10th order). Therefore, the new wavefront can be expressedas

$\begin{matrix}\begin{matrix}{{W^{\prime}\left( {u^{\prime},v^{\prime}} \right)} = {{W\left( {u,v} \right)} - {\frac{d}{2\; R^{2}}{a\left( {u,v} \right)}}}} \\{= {{\sum\limits_{i = 1}^{J}{\alpha_{i}{T_{i}\left( {u,v} \right)}}} - {\frac{d}{2\; R^{2}}{\sum\limits_{i = 1}^{J^{\prime}}{\beta_{i}{T_{i}\left( {u,v} \right)}}}}}} \\{{= {\sum\limits_{i = 1}^{J^{\prime}}{\left( {\alpha_{i} - {\frac{d}{2\; R^{2}}\beta_{i}}} \right){T_{i}\left( {u,v} \right)}}}},}\end{matrix} & ({A16})\end{matrix}$

where α_(i)=0 for i>J. With Eq. (A16), the original wavefront can beconverted from Zernike polynomials to Taylor monomials (Dai, G.-m. J.Opt. Soc. Am. A 23:1657-1666 (2006)), and the wavefront is propagatedusing Eq. (A16), then it can be converted back to Zernike polynomials.Table 6 shows the direction factor in terms of Zernike coefficientsb_(i) after propagation as a function of the coefficients c_(i) beforepropagation. The new wavefront W′(u′, v′) can be expressed within thenew wavefront boundary. For low order aberrations, since the newboundary becomes elliptical, the new elliptical wavefront can beconverted to a circular wavefront using the classical vertex correctionformula to a given new wavefront radius. For high order aberrations, thewavefront map can be be rescaled (Dai, G.-m. J. Opt. Soc. Am. A23:539-543 (2006)) to account for the change of the wavefront radius.

Table 5 relates β and α. Their relation appears to be nonlinear. In someembodiments, Taylor coefficients β_(i) of the direction factor expressedas Taylor coefficients of the original wavefront a_(i), up to the 10thorder.

TABLE 5 Symbol Expression β₁ 4α₁α₃ + 2α₂α₄ β₂ 2α₁α₄ + 4α₂α₅ β₃ 6α₁α₆ +2α₂α₇ + 4α₃ ² + α₄ ² β₄ 4α₁α₇ + 4α₂α₈ + 4α₃α₄ + 4α₄α₅ β₅ 2α₁α₈ + 6α₂α₉ +α₄ ² + 4α₅ ² β₆ 8α₁α₁₀ + 2α₂α₁₁ + 12α₃α₆ + 2α₄α₇ β₇ 6α₁α₁₁ + 4α₂α₁₂ +8α₃α₇ + 6α₄α₆ + 4α₄α₈ + 4α₅α₇ β₈ 4α₁α₁₂ + 6α₂α₁₃ + 4α₃α₈ + 4α₄α₇ +6α₄α₉ + 8α₅α₈ β₉ 2α₁α₁₃ + 8α₂α₁₄ + 2α₄α₈ + 12α₅α₉ β₁₀ 10α₁α₁₅ + 2α₂α₁₆ +16α₃α₁₀ + 2α₄α₁₁ + 9α₆ ² + α₇ ² β₁₁ 8α₁α₁₆ + 4α₂α₁₇ + 12α₃α₁₁ + 8α₄α₁₀ +4α₄α₁₂ + 4α₅α₁₁ + 12α₆α₇ + 4α₇α₈ β₁₂ 6α₁α₁₇ + 6α₂α₁₈ + 8α₃α₁₂ + 6α₄α₁₁ +6α₄α₁₃ + 8α₅α₁₂ + 6α₆α₈ + 4α₇ ² + 6α₇α₉ + 4α₈ ² β₁₃ 4α₁α₁₈ + 8α₂α₁₉ +4α₃α₁₃ + 4α₄α₁₂ + 8α₄α₁₄ + 12α₅α₁₃ + 4α₇α₈ + 12α₈α₉ β₁₄ 2α₁α₁₉ +10α₂α₂₀ + 2α₄α₁₃ + 16α₅α₁₄ + α₈ ² + 9α₉ ² β₁₅ 12α₁α₂₁ + 2α₂α₂₂ +20α₃α₁₅ + 2α₄α₁₆ + 24α₆α₁₀ + 2α₇α₁₁ β₁₆ 10α₁α₂₂ + 4α₂α₂₃ + 16α₃α₁₆ +10α₄α₁₅ + 4α₄α₁₇ + 4α₅α₁₆ + 18α₆α₁₁ + 16α₇α₁₀ + 4α₇α₁₂ + 4α₈α₁₁ β₁₇8α₁α₂₃ + 6α₂α₂₄ + 12α₃α₁₇ + 8α₄α₁₆ + 6α₄α₁₈ + 8α₅α₁₇ + 12α₆α₁₂ +12α₇α₁₁ + 6α₇α₁₃ + 8α₈α₁₀ + 8α₈α₁₂ + 6α₉α₁₁ β₁₈ 6α₁α₂₄ + 8α₂α₂₅ +8α₃α₁₈ + 6α₄α₁₇ + 8α₄α₁₉ + 12α₅α₁₈ + 6α₆α₁₃ + 8α₇α₁₂ + 8α₇α₁₄ + 6α₈α₁₁ +12α₈α₁₃ + 6α₉α₁₂ β₁₉ 4α₁α₂₅ + 10α₂α₂₆ + 4α₃α₁₉ + 4α₄α₁₈ + 10α₄α₂₀ +16α₅α₁₉ + 4α₇α₁₃ + 4α₈α₁₂ + 16α₈α₁₄ + 18α₉α₁₃ β₂₀ 2α₁α₂₆ + 12α₂α₂₇ +2α₄α₁₉ + 20α₅α₂₀ + 2α₈α₁₃ + 24α₉α₁₄ β₂₁ 14α₁α₂₈ + 24α₃α₂₁ + 2α₄α₂₂ +30α₆α₁₅ + 2α₇α₁₆ + 16α₁₀ ² + α₁₁ ² β₂₂ 20α₃α₂₂ + 12α₄α₂₁ + 4α₄α₂₃ +4α₅α₂₂ + 24α₆α₁₆ + 20α₇α₁₅ + 4α₇α₁₇ + 4α₈α₁₆ + 24α₁₀α₁₁ + 4α₁₁α₁₂ β₂₃16α₃α₂₃ + 10α₄α₂₂ + 6α₄α₂₄ + 8α₅α₂₃ + 18α₆α₁₇ + 16α₇α₁₆ + 6α₇α₁₈ +10α₈α₁₅ + 8α₈α₁₇ + 6α₉α₁₆ + 16α₁₀α₁₂ + 9α₁₁ ² + 6α₁₁α₁₃ + 4α₁₂ ² β₂₄12α₃α₂₄ + 8α₄α₂₃ + 8α₄α₂₅ + 12α₅α₂₄ + 12α₆α₁₈ + 12α₇α₁₇ + 8α₇α₁₉ +8α₈α₁₆ + 12α₈α₁₈ + 12α₉α₁₇ + 8α₁₀α₁₃ + 12α₁₁α₁₂ + 8α₁₁α₁₄ + 12α₁₂α₁₃ β₂₅8α₃α₂₅ + 6α₄α₂₄ + 10α₄α₂₆ + 16α₅α₂₅ + 6α₆α₁₉ + 8α₇α₁₈ + 10α₇α₂₀ +6α₈α₁₇ + 16α₈α₁₉ + 18α₉α₁₈ + 6α₁₁α₁₃ + 4α₁₂ ² + 16α₁₂α₁₄ + 9α₁₃ ² β₂₆4α₃α₂₆ + 4α₄α₂₅ + 12α₄α₂₇ + 20α₅α₂₆ + 4α₇α₁₉ + 4α₈α₁₈ + 20α₈α₂₀ +24α₉α₁₉ + 4α₁₂α₁₃ + 24α₁₃α₁₄ β₂₇ 2α₄α₂₆ + 24α₅α₂₇ + 2α₈α₁₉ + 30α₉α₂₀ +α₁₃ ² + 16α₁₄ ² β₂₈ 28α₃α₂₈ + 36α₆α₂₁ + 2α₇α₂₂ + 40α₁₀α₁₅ + 2α₁₁α₁₆ β₂₉14α₄α₂₈ + 30α₆α₂₂ + 24α₇α₂₁ + 4α₇α₂₃ + 4α₈α₂₂ + 32α₁₀α₁₆ + 30α₁₁α₁₅ +4α₁₁α₁₇ + 4α₁₂α₁₆ β₃₀ 24α₆α₂₃ + 20α₇α₂₂ + 6α₇α₂₄ + 12α₈α₂₁ + 8α₈α₂₃ +6α₉α₂₂ + 24α₁₀α₁₇ + 24α₁₁α₁₆ + 6α₁₁α₁₈ 20α₁₂α₁₅ + 8α₁₂α₁₇ + 6α₁₃α₁₆ β₃₁18α₆α₂₄ + 16α₇α₂₃ + 8α₇α₂₅ + 10α₈α₂₂ + 12α₈α₂₄ + 12α₉α₂₃ + 16α₁₀α₁₈ +18α₁₁α₁₇ + 8α₁₁α₁₉ 16α₁₂α₁₆ + 12α₁₂α₁₈ + 10α₁₃α₁₅ + 12α₁₃α₁₇ + 8α₁₄α₁₆β₃₂ 12α₆α₂₅ + 12α₇α₂₄ + 10α₇α₂₆ + 8α₈α₂₃ + 16α₈α₂₅ + 18α₉α₂₄ + 8α₁₀α₁₉ +12α₁₁α₁₈ + 10α₁₁α₂₀ 12α₁₂α₁₇ + 16α₁₂α₁₉ + 8α₁₃α₁₆ + 18α₁₃α₁₈ + 16α₁₄α₁₇β₃₃ 6α₆α₂₆ + 8α₇α₂₅ + 12α₇α₂₇ + 6α₈α₂₄ + 20α₈α₂₆ + 24α₉α₂₅ + 6α₁₁α₁₉ +8α₁₂α₁₈ + 20α₁₂α₂₀ 6α₁₃α₁₇ + 24α₁₃α₁₉ + 24α₁₄α₁₈ β₃₄ 4α₇α₂₆ + 4α₈α₂₅ +24α₈α₂₇ + 30α₉α₂₆ + 4α₁₂α₁₉ + 4α₁₃α₁₈ + 30α₁₃α₂₀ + 32α₁₄α₁₉ β₃₅ 2α₈α₂₆ +36α₉α₂₇ + 2α₁₃α₁₉ + 40α₁₄α₂₀ β₃₆ 42α₆α₂₈ + 48α₁₀α₂₁ + 2α₁₁α₂₂ + 25α₁₅² + α₁₆ ² β₃₇ 28α₇α₂₈ + 40α₁₀α₂₂ + 36α₁₁α₂₁ + 4α₁₁α₂₃ + 4α₁₂α₂₂ +40α₁₅α₁₆ + 4α₁₆α₁₇ β₃₈ 14α₈α₂₈ + 32α₁₀α₂₃ + 30α₁₁α₂₂ + 6α₁₁α₂₄ +24α₁₂α₂₁ + 8α₁₂α₂₃ + 6α₁₃α₂₂ + 30α₁₅α₁₇ + 16α₁₆ ² + 6α₁₆α₁₈ + 4α₁₇ ² β₃₉24α₁₀α₂₄ + 24α₁₁α₂₃ + 8α₁₁α₂₅ + 20α₁₂α₂₂ + 12α₁₂α₂₄ + 12α₁₃α₂₁ +12α₁₃α₂₃ + 8α₁₄α₂₂ + 20α₁₅α₁₈ + 24α₁₆α₁₇ + 8α₁₆α₁₉ + 12α₁₇α₁₈ β₄₀16α₁₀α₂₅ + 18α₁₁α₂₄ + 10α₁₁α₂₆ + 16α₁₂α₂₃ + 16α₁₂α₂₅ + 10α₁₃α₂₂ +18α₁₃α₂₄ + 16α₁₄α₂₃ + 10α₁₅α₁₉ + 16α₁₆α₁₈ + 10α₁₆α₂₀ + 9α₁₇ ² +16α₁₇α₁₉ + 9α₁₈ ² β₄₁ 8α₁₀α₂₆ + 12α₁₁α₂₅ + 12α₁₁α₂₇ + 12α₁₂α₂₄ +20α₁₂α₂₆ + 8α₁₃α₂₃ + 24α₁₃α₂₅ + 24α₁₄α₂₄ + 8α₁₆α₁₉ + 12α₁₇α₁₈ +20α₁₇α₂₀ + 24α₁₈α₁₉ β₄₂ 6α₁₁α₂₆ + 8α₁₂α₂₅ + 24α₁₂α₂₇ + 6α₁₃α₂₄ +30α₁₃α₂₆ + 32α₁₄α₂₅ + 6α₁₇α₁₉ + 4α₁₈ ² + 30α₁₈α₂₀ + 16α₁₉ ² β₄₃4α₁₂α₂₆ + 4α₁₃α₂₅ + 36α₁₃α₂₇ + 40α₁₄α₂₆ + 4α₁₈α₁₉ + 40α₁₉α₂₀ β₄₄2α₁₃α₂₆ + 48α₁₄α₂₇ + α₁₉ ² + 25α₂₀ ² β₄₅ 56α₁₀α₂₈ + 60α₁₅α₂₁ + 2α₁₆α₂₂β₄₆ 42α₁₁α₂₈ + 50α₁₅α₂₂ + 48α₁₆α₂₁ + 4α₁₆α₂₃ + 4α₁₇α₂₂ β₄₇ 28α₁₂α₂₈ +40α₁₅α₂₃ + 40α₁₆α₂₂ + 6α₁₆α₂₄ + 36α₁₇α₂₁ + 8α₁₇α₂₃ + 6α₁₈α₂₂ β₄₈14α₁₃α₂₈ + 30α₁₅α₂₄ + 32α₁₆α₂₃ + 8α₁₆α₂₅ + 30α₁₇α₂₂ + 12α₁₇α₂₄ +24α₁₈α₂₁ + 12α₁₈α₂₃ + 8α₁₉α₂₂ β₄₉ 20α₁₅α₂₅ + 24α₁₆α₂₄ + 10α₁₆α₂₆ +24α₁₇α₂₃ + 16α₁₇α₂₅ + 20α₁₈α₂₂ + 18α₁₈α₂₄ + 12α₁₉α₂₁ + 16α₁₉α₂₃ +10α₂₀α₂₂ β₅₀ 10α₁₅α₂₆ + 16α₁₆α₂₅ + 12α₁₆α₂₇ + 18α₁₇α₂₄ + 20α₁₇α₂₆ +16α₁₈α₂₃ + 24α₁₈α₂₅ + 10α₁₉α₂₂ + 24α₁₉α₂₄ + 20α₂₀α₂₃ β₅₁ 8α₁₆α₂₆ +12α₁₇α₂₅ + 24α₁₇α₂₇ + 12α₁₈α₂₄ + 30α₁₈α₂₆ + 8α₁₉α₂₃ + 32α₁₉α₂₅ +30α₂₀α₂₄ β₅₂ 6α₁₇α₂₆ + 8α₁₈α₂₅ + 36α₁₈α₂₇ + 6α₁₉α₂₄ + 40α₁₉α₂₆ +40α₂₀α₂₅ β₅₃ 4α₁₈α₂₆ + 4α₁₉α₂₅ + 48α₁₉α₂₇ + 50α₂₀α₂₆ β₅₄ 2α₁₉α₂₆ +60α₂₀α₂₇ β₅₅ 70α₁₅α₂₈ + 36α₂₁ ² + α₂₂ ² β₅₆ 56α₁₆α₂₈ + 60α₂₁α₂₂ +4α₂₂α₂₃ β₅₇ 42α₁₇α₂₈ + 48α₂₁α₂₃ + 25α₂₂ ² + 6α₂₂α₂₄ + 4α₂₃ ² β₅₈28α₁₈α₂₈ + 36α₂₁α₂₄ + 40α₂₂α₂₃ + 8α₂₂α₂₅ + 12α₂₃α₂₄ β₅₉ 14α₁₉α₂₈ +24α₂₁α₂₅ + 30α₂₂α₂₄ + 10α₂₂α₂₆ + 16α₂₃ ² + 16α₂₃α₂₅ + 9α₂₄ ² β₆₀12α₂₁α₂₆ + 20α₂₂α₂₅ + 12α₂₂α₂₇ + 24α₂₃α₂₄ + 20α₂₃α₂₆ + 24α₂₄α₂₅ β₆₁10α₂₂α₂₆ + 16α₂₃α₂₅ + 24α₂₃α₂₇ + 9α₂₄ ² + 30α₂₄α₂₆ + 16α₂₅ ² β₆₂8α₂₃α₂₆ + 12α₂₄α₂₅ + 36α₂₄α₂₇ + 40α₂₅α₂₆ β₆₃ 6α₂₄α₂₆ + 4α₂₅ ² +48α₂₅α₂₇ + 25α₂₆ ² β₆₄ 4α₂₅α₂₆ + 60α₂₆α₂₇ β₆₅ α₂₆ ² + 36α₂₇ ²

Note, the relationship of J′ and J is as such: in terms of degree p forJ, the new degree is p′=2p−2 for J′. For example, if J=27 (6^(th)order), then J′=65 (10^(th) order). Therefore, the new wavefront can beexpressed as

$\begin{matrix}\begin{matrix}{{W^{\prime}\left( {x^{\prime},y^{\prime}} \right)} = {{W\left( {x,y} \right)} - {\frac{d}{2\; R^{2}}\left\lbrack {\left( \frac{\partial{W\left( {x,y} \right)}}{\partial x} \right)^{2} + \left( \frac{\partial{W\left( {x,y} \right)}}{\partial y} \right)^{2}} \right\rbrack}}} \\{= {{\sum\limits_{i = 1}^{J}{\alpha_{i}{T_{i}\left( {x,y} \right)}}} - {\frac{d}{2\; R^{2}}{\sum\limits_{i = 1}^{J^{\prime}}{\beta_{i}{T_{i}\left( {x,y} \right)}}}}}} \\{= {\sum\limits_{i = 1}^{J^{\prime}}{\left( {\alpha_{i} - {\frac{d}{2\; R^{2}}\beta_{i}}} \right){T_{i}\left( {x,y} \right)}}}} \\{{= {\sum\limits_{i = 1}^{J^{\prime}}{t_{i}{T_{i}\left( {x,y} \right)}}}},}\end{matrix} & (95)\end{matrix}$

where the new coefficients t_(i) is related to the original coefficientsα_(i) and the derived coefficients β_(i) as

$\begin{matrix}{t_{i} = {\alpha_{i} - {\frac{d}{2\; R^{2}}{\beta_{i}.}}}} & (96)\end{matrix}$

Once the derived coefficients β_(i) are calculated, the new coefficientst_(i) are known. Because Taylor coefficients can be converted to Zernikecoefficients as

$\begin{matrix}{c_{i} = {\sum\limits_{j = 1}^{J^{\prime}}{C_{ij}^{t\; 2\; z}t_{j}}}} & (97)\end{matrix}$

the wavefront can be represented as Zernike polynomials as

$\begin{matrix}{{W^{\prime}\left( {\rho^{\prime},\theta^{\prime}} \right)} = {\sum\limits_{i = 1}^{J^{\prime}}{c_{i}{{Z_{i}\left( {\rho,\theta} \right)}.}}}} & (98)\end{matrix}$

Comparison of Eqs. (88) and (92), we can relate the coefficients ofα_(i) and a_(i) as

$\begin{matrix}{\alpha_{i} = {\sum\limits_{j = 1}^{J}{C_{ij}^{z\; 2\; t}a_{j}}}} & (99)\end{matrix}$

Because c_(i,)is a function of t_(i), t_(i) is a function of α_(i) andβ_(i), and β_(i) is a function of α_(i), we reason that c_(i) can becalculated from α_(i).

As shown here, Table 6 provides Zernike coefficients b_(i) of thedirection factor expressed as those in the original wavefront c_(i), upto the 6th order.

TABLE 6 Symbol Expression b₁ 4({square root over (6)}c₂c₃ + c₁(2{squareroot over (5)}c₁₂ − {square root over (10)}c₁₃ + 2{square root over(3)}c₄ − {square root over (6)}c₅) + 6{square root over (5)}c₁₃c₆ −6{square root over (5)}c₁₄c₆ + 4{square root over (3)}c₅c₆ + 14{squareroot over (10)}c₁₂c₇ − 14{square root over (5)}c₁₃c₇ + 8{square rootover (6)}c₄c₇ − 4{square root over (3)}c₅c₇ + 4{square root over(3)}c₃c₈ + {square root over (5)}c₁₁({square root over (2)}c₂ + 14c₈ −6c₉) + 6{square root over (5)}c₁₀c₉ − 4{square root over (3)}c₃c₉) b₂4({square root over (10)}c₁₃c₂ + c₁({square root over (10)}c₁₁ + {squareroot over (6)}c₃) + 2{square root over (3)}c₂c₄ + {square root over(6)}c₂c₅ + 6{square root over (5)}c₁₀c₆ + 6{square root over (5)}c₁₁c₆ +4{square root over (3)}c₃c₆ + 14{square root over (5)}c₁₁c₇ + 4{squareroot over (3)}c₃c₇ + 14{square root over (5)}c₁₃c₈ + 8{square root over(6)}c₄c₈ + 4{square root over (3)}c₅c₈ + 2{square root over (5)}c₁₂(c₂ +7{square root over (2)}c₈) + 6{square root over (5)}c₁₃c₉ + 6{squareroot over (5)}c₁₄c₉ + 4{square root over (3)}c₅c₉) b₃ 4(9{square rootover (6)}c₁₀c₁₃ + 6{square root over (5)}c₁₂c₃ − 3{square root over(10)}c₁₄c₃ + 4{square root over (3)}c₃c₄ + 3c₁₁(14{square root over(3)}c₁₂ − 3{square root over (6)}c₁₄ + 4{square root over (5)}c₄) +3{square root over (10)}c₁₀c₅ + 2{square root over (3)}c₂c₆ + 2{squareroot over (3)}c₂c₇ + 2{square root over (3)}c₁c₈ + 5{square root over(6)}c₆c₈ + 10{square root over (6)}c₇c₈ − 2{square root over (3)}c₁c₉ −5{square root over (6)}c₇c₉) b₄ 4(6{square root over (3)}c₁₀ ² +15{square root over (3)}c₁₁ ² + 18{square root over (3)}c₁₂ ² +15{square root over (3)}c₁₃ ² + 6{square root over (3)}c₁₄ ² + 6{squareroot over (5)}c₁₁c₃ + {square root over (3)}c₃ ² + 12{square root over(5)}c₁₂c₄ + 2{square root over (3)}c₄ ² + 6{square root over (5)}c₁₃c₅ +{square root over (3)}c₅ ² + 3{square root over (3)}c₆ ² + 2{square rootover (6)}c₁c₇ + 7{square root over (3)}c₇ ² + 2{square root over(6)}c₂c₈ + 7{square root over (3)}c₈ ² + 3{square root over (3)}c₉ ²) b₅4(9{square root over (6)}c₁₀c₁₁ + 42{square root over (3)}c₁₂c₁₃ +9{square root over (6)}c₁₃c₁₄ + 3{square root over (10)}c₁₀c₃ +12{square root over (5)}c₁₃c₄ + 6{square root over (5)}c₁₂c₅ + 3{squareroot over (10)}c₁₄c₅ + 4{square root over (3)}c₄c₅ + 2{square root over(3)}c₁c₆ − 2{square root over (3)}c₁c₇ + 5{square root over (6)}c₆c₇ −5{square root over (6)}c₇ ² + 2{square root over (3)}c₂c₈ + 5{squareroot over (6)}c₈ ² + 2{square root over (3)}c₂c₉ + 5{square root over(6)}c₈c₉) b₆ 8{square root over (5)}c₁(c₁₃ − c₁₄) + 4(10{square rootover (5)}c₁₁c₂ + 54{square root over (5)}c₁₂c₆ + 30{square root over(3)}c₄c₆ + 55{square root over (10)}c₁₃c₇ − 28{square root over(10)}c₁₄c₇ + 15{square root over (6)}c₅c₇ + 55{square root over(10)}c₁₁c₈ + 15{square root over (6)}c₃c₈ + 2{square root over(5)}c₁₀(5c₂ + 14{square root over (2)}c₈))/5 b₇ 4(10{square root over(5)}c₁({square root over (2)}c₁₂ − c₁₃) + 21{square root over (10)}c₁₃c₆− 12{square root over (10)}c₁₄c₆ + 5{square root over (6)}c₅c₆ +86{square root over (5)}c₁₂c₇ − 34{square root over (10)}c₁₃c₇ +30{square root over (3)}c₄c₇ − 10{square root over (6)}c₅c₇ + 10{squareroot over (6)}c₃c₈ + {square root over (5)}c₁₁(10c₂ + {square root over(2)}(34c₈ − 21c₉)) + 12{square root over (10)}c₁₀c₉ − 5{square root over(6)}₃c₉)/5 b₈ 4(10{square root over (5)}c₁c₁₁ + 10{square root over(5)}c₁₃c₂ + 12{square root over (10)}c₁₀c₆ + 21{square root over(10)}c₁₁c₆ + 5{square root over (6)}c₃c₆ + 34{square root over(10)}c₁₁c₇ + 10{square root over (6)}c₃c₇ + 34{square root over(10)}c₁₃c₈ + 30{square root over (3)}c₄c₈ + 10{square root over(6)}c₅c₈ + 2{square root over (5)}c₁₂(5{square root over (2)}c₂ +43c₈) + 21{square root over (10)}c₁₃c₉ + 12{square root over(10)}c₁₄c₉ + 5{square root over (6)}c₅c₉)/5 b₉ 8{square root over(5)}c₁(c₁₀ − c₁₁) + 4(10{square root over (5)}c₁₄c₂ + 28{square rootover (10)}c₁₀c₇ − 55{square root over (10)}c₁₁c₇ − 15{square root over(6)}c₃c₇ + 28{square root over (10)}c₁₄c₈ + 15{square root over(6)}c₅c₈ + 5{square root over (5)}c₁₃(2c₂ + 11{square root over(2)}c₈) + 54{square root over (5)}c₁₂c₉ + 30{square root over(3)}c₄c₉)/5 b₁₀ 8(10{square root over (6)}₁₃c₃ + 20c₁₀(2{square rootover (5)}c₁₂ + {square root over (3)}c₄) + 10c₁₁(4{square root over(10)}c₁₃ + {square root over (6)}c₅) + 9{square root over (10)}c₆c₈ +9{square root over (10)}c₇c₉)/5 b₁₁ 4(30{square root over (3)}c₁₂c₃ −5{square root over (6)}c₁₄c₃ + 5c₁₁(22{square root over (5)}c₁₂ −5{square root over (10)}c₁₄ + 8{square root over (3)}c₄) + 5c₁₀(5{squareroot over (10)}c₁₃ + {square root over (6)}c₅) + 9{square root over(10)}c₆c₈ + 18{square root over (10)}c₇c₈ − 9{square root over(10)}c₇c₉)/5 b₁₂ 8{square root over (5)}c₁₀ ² + 32{square root over(5)}c₁₁ ² + 48{square root over (5)}c₁₂ ² + 32{square root over (5)}c₁₃² + 8{square root over (5)}c₁₄ ² + 16{square root over (3)}c₁₁c₃ +32{square root over (3)}c₁₂c₄ + 16{square root over (3)}c₁₃c₅ + (12c₆²)/{square root over (5)}+ 12{square root over (5)}c₇ ² + 12{square rootover (5)}c₈ ² + (12c₉ ²)/{square root over (5)} b₁₃ 4(110{square rootover (5)}c₁₂c₁₃ + 25{square root over (10)}c₁₃c₁₄ + 5c₁₀(5{square rootover (10)}c₁₁ + {square root over (6)}c₃) + 40{square root over(3)}c₁₃c₄ + 30{square root over (3)}₁₂c₅ + 5{square root over(6)}c₁₄c₅ + 9{square root over (10)}c₆c₇ − 9{square root over (10)}c₇² + 9{square root over (10)}c₈ ² + 9{square root over (10)}c₈c₉)/5 b₁₄−32{square root over (10)}c₁₁ ² + 32{square root over (10)}c₁₃ ² +64{square root over (5)}c₁₂c₁₄ − 16{square root over (6)}c₁₁c₃ +32{square root over (3)}c₁₄c₄ + 16{square root over (6)}c₁₃c₅ −72{square root over (2/5)}c₆c₇ + 72{square root over (2/5)}c₈c₉ b₁₅16{square root over (15)}(c₁₃c₆ + c₁₄c₇ + c₁₀c₈ + c₁₁c₉) b₁₆ 16{squareroot over (3/5)}(3{square root over (2)}c₁₂c₆ + 5c₁₃c₇ − 2c₁₄c₇ +2c₁₀c₈ + 5c₁₁c₈) b₁₇ 8{square root over (3/5)}(3c₁₃c₆ − c₁₄c₆ + 9{squareroot over (2)}c₁₂c₇ − 7c₁₃c₇ + 7c₁₁c₈ + c₁₀c₉ − 3c₁₁c₉) b₁₈ 8{squareroot over (3/5)}(c₁₀c₆ + 3c₁₁c₆ + 7c₁₁c₇ + 9{square root over(2)}c₁₂c₈ + 7c₁₃c₈ + 3c₁₃c₉ + c₁₄c₉) b₁₉ 16{square root over(3/5)}(2c₁₀c₇ − 5c₁₁c₇ + 5c₁₃c₈ + 2c₁₄c₈ + 3{square root over (2)}c₁₂c₉)b₂₀ −16{square root over (15)}(c₁₁c₆ + c₁₀c₇ − c₁₄c₈ − c₁₃c₉) b₂₁160{square root over (2/7)}(c₁₀c₁₃ + {square root over (2)}c₁₁c₁₄) b₂₂160(c₁₀c₁₂ + {square root over (2)}c₁₁c₁₃)/{square root over (7)} b₂₃32(8c₁₁c₁₂ + {square root over (2)}c₁₀c₁₃ − {square root over(2)}c₁₁c₁₄)/{square root over (7)} b₂₄ 8(c₁₀ ² + 10c₁₁ ² + 18c₁₂ ² +10c₁₃ ² + c₁₄ ²)/{square root over (7)} b₂₅ 32({square root over(2)}c₁₀c₁₁ + 8c₁₂c₁₃ + {square root over (2)}c₁₃c₁₄)/{square root over(7)} b₂₆ −80({square root over (2)}c₁₁ ² − {square root over (2)}c₁₃ ² −2c₁₂c₁₄)/{square root over (7)} b₂₇ −160{square root over (2/7)}(c₁₀c₁₁− c₁₃c₁₄)

Table 7 provides a comparison of Zernike Coefficients, in μm, after alow order wavefront propagated by 3.5 mm or 12.5 mm using the vertexcorrection method and Zemax®. Both the original and propagatedwavefronts are represented over a 6 mm pupil.

TABLE 7 Output (d = 3.5 mm) Out (d = 12.5 mm) Cost Term Input VertexZemax ® Diff (%) Vertex Zemax ® Diff (%) One c₂ ⁻² 0.3140 0.3086 0.30850.03% 0.2953 0.2944 0.30% c₂ ⁰ −3.2480 −3.2187 −3.2181 0.02% −3.1459−3.1420 0.12% c₂ ² −0.8630 −0.8481 −0.8503 0.26% −0.8115 −0.8089 0.32%Two c₂ ⁻² −2.2910 −2.2951 −2.2950 0.00% −2.3060 −2.3053 0.03% c₂ ⁰0.3250 0.3324 0.3323 0.03% 0.3516 0.3514 0.06% c₂ ² 0.1600 0.1603 0.16020.05% 0.1610 0.1610 0.00%

Table 8 provides a comparison of Zernike coefficients, in μm, after asingle term wavefront ropagated by 12.5 mm using the analytical methodand Zemax®. Both the original and propagated wavefronts are representedover a 6 mm pupil.

TABLE 8 Term Analytical Zemax ® c₄ ⁰ = 0.5 μm c₂ ⁰ −0.217 −0.217 c₄ ⁰0.4814 0.4813 c₂ ⁰ −0.0094 −0.0096 c₃ ¹ = 0.5 μm c₂ ⁰ −0.0084 −0.0180 c₂² −0.0085 −0.0153 c₃ ¹ 0.5000 0.5000 c₄ ⁰ −0.0047 −0.0047 c₄ ² −0.0040−0.0040 c₃ ³ = 0.5 μm c₂ ⁰ −0.0036 −0.0036 c₃ ³ 0.5000 0.5000 c₄ ⁰−0.0009 −0.0009

Table 9 provides a comparison of Zernike coefficients, in μm, after arandom wavefront propagated by 12.5 mm using the analytical method andZemax®. Both the original and propagated wavefronts are represented overa 6 mm pupil.

TABLE 9 Output Term Input Analytical Zemax ® Diff c₂ ⁻² −0.0646 −0.0757−0.0795 −0.0038 c₂ ⁰ 3.3435 3.4567 3.4588 0.0021 c₂ ² −0.0217 −0.0230−0.0240 −0.0010 c₃ ⁻³ −0.2103 −0.2222 −0.2215 0.0007 c₃ ⁻¹ 0.0810 0.08550.0786 −0.0069 c₃ ¹ 0.0590 0.0659 0.0763 0.0104 c₃ ³ 0.0047 0.0018−0.0018 −0.0036 c₄ ⁻⁴ 0.0451 0.0583 0.0541 −0.0042 c₄ ⁻² −0.0538 −0.0565−0.0569 −0.0042 c₄ ⁰ 0.0705 0.0732 0.0918 0.0186 c₄ ² 0.1110 0.11020.1236 0.0134 c₄ ⁴ −0.0477 −0.0535 −0.0631 −0.0096 c₅ ⁻⁵ 0.0169 0.02110.0206 −0.0005 c₅ ⁻³ 0.0276 0.0263 0.0315 0.0052 c₅ ⁻¹ −0.0296 −0.0292−0.0344 −0.0052 c₅ ¹ 0.0294 0.0315 0.0358 0.0043 c₅ ³ −0.0217 −0.0231−0.0257 −0.0026 c₅ ⁵ −0.0079 −0.0067 −0.0102 −0.0035 c₆ ⁻⁶ 0.0575 0.06010.0706 0.0105 c₆ ⁻⁴ 0.0080 0.0111 0.0102 −0.0009 c₆ ⁻² 0.0122 0.01230.0141 0.0018 c₆ ⁰ 0.0280 0.0286 0.0340 0.0054 c₆ ² −0.0093 −0.0109−0.0101 0.0008 c₆ ⁴ −0.0179 −0.0190 −0.0233 −0.0043 c₆ ⁶ −0.0332 −0.0360−0.0405 −0.0045 r.m.s 11.2716 12.0531 12.0785 0.0331 ho rms −0.08790.0980 0.1083 0.0320

FIG. 13 shows a flow chart for the process for an original wavefrontpropagated to a certain distance to become a propagated new wavefront,using an analytical approach. The original wavefront is expanded intoZernike polynomials using Eq. (88). The first conversion into Taylormonomials uses Eq. (92), with the conversion matrix given by Eq. (A).

$\begin{matrix}{{C_{z\; 2\; t} = {\frac{\left( {- 1} \right)^{{({n - q})}/2}{{{{\sqrt{n + 1}\left\lbrack {\left( {n + p} \right)/2} \right\rbrack}!}\lbrack m\rbrack}!}}{{{\left\lbrack {\left( {p + {m}} \right)/2} \right\rbrack!}\left\lbrack {\left( {n - p} \right)/2} \right\rbrack}!} \times {\sum\limits_{t = 0}^{t_{0}}{\sum\limits_{t^{\prime} = 0}^{{({p - {m}})}/2}\frac{\left( {- 1} \right)^{t}\sqrt{2 - \delta_{m\; 0}}{f\left( {m,t} \right)}}{{{\left( t^{\prime} \right)!}\left\lbrack {{\left( {p - {m}} \right)/2} - t^{\prime}} \right\rbrack}!}}}}},} & (A)\end{matrix}$

The propagation is performed with Eq. (96) using the propagationdistance d, and the new wavefront diameter is obtained with Eq. (87).The second conversion, which converts Taylor coefficients to Zernikecoefficients, is done with Eq. (97) using the conversion matrix fromgiven by Eq. (B). Finally, the new Zernike coefficients and newwavefront diameter are obtained, and the combination gives the newpropagated wavefront.

$\begin{matrix}\begin{matrix}{C_{t\; 2\; z} = {\frac{1}{\pi}{\int_{0}^{1}{{\int_{0}}^{2\; \pi}{{T_{p}^{q}\left( {\rho,\theta} \right)}{Z_{n}^{m}\left( {\rho,\theta} \right)}\ {\rho}}}}}} \\{= {\frac{{\left( {p - q} \right)!}{q!}}{2^{p}}{\sum\limits_{t - 0}^{q}{\frac{1}{{t!}{\left( {q - t} \right)!}}{\sum\limits_{t^{\prime} = 0}^{p - q}{\frac{g\left( {p,q,m,t,t^{\prime}} \right)}{{\left( t^{\prime} \right)!}{\left( {p - q - t^{\prime}} \right)!}} \times}}}}}} \\{{\sum\limits_{s = 0}^{{({n - {m}})}/2}{\frac{\left( {- 1} \right)^{s}\sqrt{n + 1}{\left( {n - s} \right)!}}{{{{{s!}\left\lbrack {{\left( {n + m} \right)/2} - s} \right\rbrack}!}\left\lbrack {{\left( {n - m} \right)/2} - s} \right\rbrack}!} \times}}} \\{{\frac{1}{n + p - {2\; s} + 2},}}\end{matrix} & (B)\end{matrix}$

G. General Numerical Representation

For a pure numerical approach, the original wavefront may be in ananalytical form, such as represented as Zernike polynomials or otherbasis functions. It can be sampled in discrete form, i.e., to obtainvalues in a 2-D space over the wavefront diameter. If the originalwavefront is already discrete, such as from a Fourier reconstruction, wecan keep the same discrete data form or we can do a resampling.Calculation of the x and y derivatives can be done by calculating thedifference of two neighboring points in either x or y direction, dividedby the sampling rate, or the distance between the two neighboringpoints. To obtain the derivative to r, the following simple formula canbe used:

$\begin{matrix}{\frac{\partial{W\left( {r,\theta} \right)}}{\partial r} = \sqrt{\left\lbrack \frac{\partial{W\left( {r,\theta} \right)}}{\partial x} \right\rbrack^{2} + \left\lbrack \frac{\partial{W\left( {r,\theta} \right)}}{\partial y} \right\rbrack^{2}}} & (100)\end{matrix}$

Once the derivative to r is calculated, the new wavefront diameter canbe obtained using Eq. (87). The new discrete wavefront can then beobtained with Eq. (38). FIG. 14 illustrates a flow chart for the processof wavefront propagation using numerical approach.

VI. Verification

For low order aberrations, it may be helpful to consider FIG. 5. Thesimple geometry shows

$\begin{matrix}{{\frac{R}{R^{\prime}} = \frac{f}{f + d}}{or}} & ({A30}) \\{R^{\prime} = {\left( {1 + \frac{d}{f}} \right){R.}}} & ({A31})\end{matrix}$

Using a plus cylinder notation, we first consider the meridian of theminimum power S. Because the focal length is related to the refractivepower by f=1/S, we obtain the semiminor axis of the propagated wavefrontas

$\begin{matrix}\begin{matrix}{R_{\min} = {\left( {1 + {dS}} \right)R}} \\{{= {R\left\{ {1 - {\frac{d}{R^{2}}\left\lbrack {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}} \right\}}},}\end{matrix} & ({A32})\end{matrix}$

as the sphere power S is given by Eq. (A56a). Similarly, for themeridian of the maximum power S+C, we have f=1/(S+C). Therefore thesemimajor axis of the propagated wavefront is

$\begin{matrix}\begin{matrix}{R_{\max} = \left\lbrack {1 + {{d\left( {S + C} \right\rbrack}R}} \right.} \\{{= {R\left\{ {1 - {\frac{d}{R^{2}}\left\lbrack {{4\sqrt{3}c_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}} \right\}}},}\end{matrix} & ({A33})\end{matrix}$

as the refractive power S and C are given by Eq. (A56). Equations (A32)and (A33) are equivalent to Eqs. (A64) and (A67), respectively.Therefore, the boundary change of the low order wavefront propagation isverified. As for the magnitude of the propagated wavefront, a detailedproof is given in Section VIII below.

For the verification of the high order aberrations, it is possible touse Zemax® software (Zemax Development Corporation, Bellevue, Wash.) asa ray tracing tool for the comparison purpose. It is possible to use thefree space propagation of a wavefront represented with Zernikepolynomials. A Hartmann-Shack wavefront sensor was attached in theZemax® model, but the calculation of the propagation of the wavefrontwas not affected by the sensor. For Zernike representation of thewavefronts from the Zemax® model, it is possible to use a 512×512wavefront sampling to reduce the fitting error. A proper Zernikecoefficient conversion can be performed as Zemax® uses the Nolls'snotation (Noll, R. J. J. Opt. Soc. Am. 66:203-211 (1976)) and the ANSInotation (American National Standard Institute, Methods for reportingoptical aberrations of eyes, ANSI Z80.28-2004 (Optical LaboratoriesAssociation, 2004), Annex B, pp. 1928) was used in embodimentsencompassed herein. In some cases, there is a slight difference for thepropagation of low order aberrations between the classical vertexcorrection formula and Zemax®, as can be seen from Table 7.

For high order aberrations, it is possible to use a few single Zernikemode aberrations as shown in Table 8. This shows these examples ofaberrations as measured on the exit pupil plane and as represented on apropagated vertex plane. A wavefront diameter of 6 mm is assumed beforeand after the propagation. Although the wavefront boundary changes afterthe propagation, the Zernike coefficients are properly scaled (Dai,G.-m. J. Opt. Soc. Am. A 23:539-543 (2006)) to the original wavefrontradius. The results are compared to those obtained with the analyticalexpressions developed in the previous section. For each of the singleZernike mode aberrations, both approaches give nearly exact results,except for low orders for the propagation of coma (Z₃ ¹). Thisdiscrepancy may be attributed to the approximation of the ellipticalpupil to a circular one.

As a further verification, it is possible to use an ocular wavefrontfrom a real eye that consists of all the first 27 Zernike coefficients,as shown in Table 9. The results obtained with the analytical approachand those obtained with Zemax® are again comparable. The differencesbetween the results using the analytical approach and Zemax® are alsoshown. Also shown are the root mean square (RMS) and high order RMSvalues. Similar to the previous examples, the wavefront radius beforeand after propagation is 6 mm. Both approaches give similar results. Thereason that the two sets of results are not identical can be attributedto approximations during the theoretical development in the previoussection as well as the numerical fitting error in Zemax® software. Inparticular, for the propagation of coma, the approximation of anelliptical pupil to a circular pupil affects the values of the induceddefocus and astigmatism. Without these approximations and numericalerror, the results can be expected to be very close, if not identical.

Embodiments encompass verification approaches that involve geometricoptics. It can be helpful to verify the wavefront-approach using aray-tracing (geometric optics) software, such as Zemax®. Beforeimplementing the software, it can be helpful to validate the Zemax®software using the sphere only and sphere and cylinder propagation.

Some validation embodiments of the present invention involving Zemax®include using a general eye model and creating a surface using Zernikepolynomials to include sphere or sphere and cylinder, i.e., Z3, Z4, andZ5. The model can be validated by adjusting a few propagation distancesto see if the ray distribution on the imaging plane is the same if thepropagated surface, in terms or Zernike polynomials, is optimized. FIG.15 shows a flow chart of a wavefront propagation using Zemax® modelingfor such a process. The result is to compare to that predicted by theclassical formula, Eqs. (58) and (59).

Once the first step is validated, the Zemax® model can then be used tovalidate the techniques described earlier. A random wavefront consistingall of the 27 terms in Zernike polynomials can be used and the resultcan be compared to those in the analytical subsection.

VII. Discussion

In some embodiments, when a wavefront propagates, it propagates as awhole and cannot be linearly combined, because the direction factor isnot a linear function, but a quadratic function, of Zernike polynomialsor Taylor monomials. Even so, it is useful to discuss some importantaberrations on how they propagate individually.

A. Low Order Aberrations

For a wavefront that consists of low order aberrations only, expressedwith Zernike polynomials, the direction factor can be written as

$\begin{matrix}{a = {{16\sqrt{3}c_{2}^{- 2}c_{2}^{0}Z_{2}^{- 2}} + {4{\sqrt{3}\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {16\sqrt{3}c_{2}^{0}c_{2}^{2}Z_{2}^{2}} + {{12\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}{Z_{0}^{0}.}}}} & ({A34})\end{matrix}$

After propagation, it does not induce any high order aberrations. Theboundary factor can be written as

$\begin{matrix}{{\left. {b^{2} = {{24\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack} + {48\sqrt{2}c_{2}^{0}}}} \right)\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}\cos \; 2\left( {\theta - \varphi} \right)},\mspace{79mu} {where}} & ({A35}) \\{\mspace{79mu} {\varphi = {\frac{1}{2}{{\tan^{- 1}\left( \frac{c_{2}^{- 2}}{c_{2}^{2}} \right)}.}}}} & ({A36})\end{matrix}$

Therefore, after propagation, the circular wavefront becomes elliptical,as shown in FIG. 12A.

B. Coma Aberration

For coma aberration, including Z₃ ⁻¹ and Z₃ ¹ Zernike terms, thedirection factor after propagation is

$\begin{matrix}{\left. {a = {{{56\left\lbrack {\left( c_{3}^{- 1} \right)^{2} + \left( c_{3}^{1} \right)^{2}} \right\rbrack}Z_{0}^{0}} + {40\sqrt{6}c_{3}^{- 1}c_{3}^{1}Z_{2}^{- 2}} + {28{\sqrt{3}\left\lbrack {\left( c_{3}^{- 1} \right)^{2} + \left( c_{3}^{1} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {20{\sqrt{6}\left\lbrack {\left( c_{3}^{1} \right)^{2} - \left( c_{3}^{- 1} \right)^{2}} \right\rbrack}Z_{2}^{2}} + {72\sqrt{\frac{2}{5}}c_{3}^{- 1}c_{3}^{1}Z_{4}^{- 2}} + {12\sqrt{5}\left( c_{3}^{- 1} \right)^{2}} + {\left( c_{3}^{1} \right)^{2}Z_{4}^{0}} + {36{\sqrt{\frac{2}{5}}\left\lbrack {\left( c_{3}^{1} \right)^{2} - \left( c_{3}^{- 1} \right)^{2}} \right\rbrack}Z_{4}^{2}}}} \right).} & ({A37})\end{matrix}$

After propagation, coma aberration induces defocus, astigmatism,spherical aberration and secondary astigmatism. The boundary factor canbe written as

$\begin{matrix}{{b^{2} = {{200\left\lbrack \left( c_{3}^{- 1} \right)^{2} \right\rbrack} + \left( c_{3}^{1} \right)^{2} + {{192\left\lbrack {\left( c_{3}^{- 1} \right)^{2} + \left( c_{3}^{1} \right)^{2} + \left( c_{3}^{1} \right)^{2}} \right\rbrack}\cos \; 2\left( {\theta - \varphi} \right)}}},\mspace{79mu} {where}} & ({A38}) \\{\mspace{79mu} {\varphi = {\frac{1}{2}{{\tan^{- 1}\left\lbrack \frac{2c_{3}^{- 1}c_{3}^{1}}{\left( c_{3}^{1} \right)^{2} - \left( c_{3}^{1} \right)^{2}} \right\rbrack}.}}}} & ({A39})\end{matrix}$

Therefore, propagation of a coma aberration becomes elliptical.

C. Trefoil Aberration

For trefoil aberration, including Z₃ ⁻³ and Z₃ ³ Zernike terms, thedirection factor after propagation becomes

$\begin{matrix}{a = {{{24\left\lbrack {\left( c_{3}^{- 3} \right)^{2} + \left( c_{3}^{3} \right)^{2}} \right\rbrack}Z_{0}^{0}} + {12{\sqrt{3}\left\lbrack {\left( c_{3}^{- 3} \right)^{2} + \left( c_{3}^{3} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {{\frac{12}{\sqrt{5}}\left\lbrack {\left( c_{3}^{- 3} \right)^{2} + \left( c_{3}^{3} \right)^{2}} \right\rbrack}{Z_{4}^{0}.}}}} & ({A40})\end{matrix}$

Therefore, propagation of trefoil only induces defocus and sphericalaberration. The boundary factor b is

b ²=72[(c ₃ ⁻³)²+(c ₃ ³)²].  (A41)

So after propagation, the boundary of an original trefoil still iscircular.

D. Primary Spherical Aberration

For the primary spherical aberration, Z₄ ⁰, the direction factor is

$\begin{matrix}{a = {{120\left( c_{4}^{0} \right)^{2}Z_{0}^{0}} + {72\sqrt{3}\left( c_{4}^{0} \right)^{2}Z_{2}^{0}} + {48\sqrt{5}\left( c_{4}^{0} \right)^{2}Z_{4}^{0}} + {\frac{144}{\sqrt{7}}\left( c_{4}^{0} \right)^{2}{Z_{6}^{0}.}}}} & ({A42})\end{matrix}$

Propagation of primary spherical aberration induces defocus, sphericalaberration and secondary spherical aberration. The boundary factor is

b ²=720(c ₄ ⁰)².  A(43)

Hence, propagation of spherical aberration still is circular.

E. Secondary Spherical Aberration

For the secondary spherical aberration, Z₆ ⁰ the direction factor is

$\begin{matrix}{a = {{336\left( c_{6}^{0} \right)^{2}Z_{0}^{0}} + {240\sqrt{3}\left( c_{6}^{0} \right)^{2}Z_{2}^{0}} + {192\sqrt{5}\left( c_{6}^{0} \right)^{2}Z_{4}^{0}} + \left. \quad{128\sqrt{7}\left( c_{6}^{0} \right)^{2}Z_{6}^{0}} \right)^{2} + {\quad{{\left( c_{6}^{0} \right)^{2}Z_{8}^{0}} + {\frac{400}{\sqrt{11}}\left( c_{6}^{0} \right)^{2}{Z_{10}^{0}.}}}}}} & ({A44})\end{matrix}$

Propagation of a secondary spherical aberration induces defocus,spherical aberration, secondary, tertiary, and quaternary sphericalaberration. The boundary factor b is

b ²=4032(c ₆ ⁰)².  (A45)

So propagation of secondary spherical aberration still is circular.

F. Secondary Astigmatism Aberration

For the secondary astigmatism, Z₄ ⁻² and Z₄ ², the direction factor is

$\begin{matrix}{\left. {\left. {a = {{100\left\lbrack {\left( c_{4}^{- 2} \right)^{2} + {\left( c_{4}^{2} \right)^{2}Z_{0}^{0}}} \right\rbrack} + {60{\sqrt{3}\left\lbrack {\left( c_{4}^{- 2} \right)^{2} + \left( c_{4}^{2} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {64\sqrt{10}c_{4}^{- 2}c_{4}^{2}Z_{4}^{- 4}} + {32{\sqrt{5}\left\lbrack {\left( c_{4}^{- 2} \right)^{2} + \left( c_{4}^{2} \right)^{2}} \right\rbrack}Z_{4}^{0}} + {32{\sqrt{10}\left\lbrack {\left( c_{4}^{2} \right)^{2} - \left( c_{4}^{- 2} \right)^{2}} \right\rbrack}Z_{4}^{4}} + {160\sqrt{\frac{2}{7}}c_{4}^{- 2}c_{4}^{2}Z_{6}^{- 4}}}} \right) + {{\frac{80}{\sqrt{7}}\left\lbrack {\left( c_{4}^{- 2} \right)^{2} + \left( c_{4}^{2} \right)^{2}} \right\rbrack}Z_{6}^{0}} + {80{\sqrt{\frac{2}{7}}\left\lbrack {\left( c_{4}^{2} \right)^{2} - \left( c_{4}^{- 2} \right)^{2}} \right\rbrack}Z_{6}^{4}}} \right).} & ({A46})\end{matrix}$

Therefore, the propagation of a secondary astigmatism induces defocus,primary and secondary spherical aberration, quadrafoil and secondaryquadrafoil. The boundary factor b is

$\begin{matrix}{{b^{2} = {{520\left\lbrack {\left( c_{4}^{- 2} \right)^{2} + \left( c_{4}^{2} \right)^{2}} \right\rbrack} + {{480\left\lbrack {\left( c_{4}^{2} \right)^{2} + \left( c_{4}^{2} \right)^{2}} \right\rbrack}\cos \; 4\left( {\theta - \varphi} \right)}}},{where}} & ({A47}) \\{\varphi = {\frac{1}{4}{{\tan^{- 1}\left\lbrack \frac{2c_{4}^{- 2}c_{4}^{2}}{\left( c_{4}^{2} \right)^{2} - \left( c_{4}^{- 2} \right)^{2}} \right\rbrack}.}}} & ({A48})\end{matrix}$

Therefore, propagation of secondary astigmatism becomes a fourfoldsymmetry of shape.

G. Secondary Coma Aberration

$\begin{matrix}{a = {{{204\left\lbrack {\left( c_{5}^{- 1} \right)^{2} + \left( c_{5}^{1} \right)^{2}} \right\rbrack}Z_{0}^{0}} + {168\sqrt{6}c_{5}^{- 1}c_{5}^{1}Z_{2}^{- 2}} + {136{\sqrt{3}\left\lbrack {\left( c_{5}^{- 1} \right)^{2} + \left( c_{5}^{1} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {84{\sqrt{6}\left\lbrack {\left( c_{5}^{- 1} \right)^{2} - \left( c_{5}^{- 1} \right)^{2}} \right\rbrack}Z_{2}^{2}} + {\frac{3576}{7}\sqrt{\frac{2}{5}}c_{5}^{- 1}c_{5}^{- 1}Z_{4}^{- 2}} + {\frac{696}{7}{\sqrt{5}\left\lbrack {\left( c_{5}^{- 1} \right)^{2} + \left( c_{5}^{1} \right)^{2}} \right\rbrack}Z_{4}^{0}} + {\frac{1788}{7}{\sqrt{\frac{2}{5}}\left\lbrack {\left( c_{5}^{1} \right)^{2} - \left( c_{5}^{- 1} \right)^{2}} \right\rbrack}Z_{4}^{2}} + {456\sqrt{\frac{2}{7}}c_{5}^{- 1}c_{5}^{1}Z_{6}^{- 2}} + {{\frac{408}{\sqrt{7}}\left\lbrack {\left( c_{5}^{- 1} \right)^{2} + \left( c_{5}^{1} \right)^{2}} \right\rbrack}Z_{6}^{0}} + {228{\sqrt{\frac{2}{7}}\left\lbrack {\left( c_{5}^{1} \right)^{2} - \left( c_{5}^{- 1} \right)^{2}} \right\rbrack}Z_{6}^{2}} + {\frac{600}{7}\sqrt{2}c_{5}^{- 1}c_{5}^{1}Z_{8}^{- 2}} + {{\frac{520}{7}\left\lbrack {\left( c_{5}^{- 1} \right)^{2} + \left( c_{5}^{1} \right)^{2}} \right\rbrack}Z_{8}^{0}} + {\frac{300}{7}{\sqrt{2}\left\lbrack {\left( c_{5}^{1} \right)^{2} - \left( c_{5}^{- 1} \right)^{2}} \right\rbrack}{Z_{8}^{2}.}}}} & ({A49})\end{matrix}$

Propagation of secondary coma induces defocus, primary, secondary, andtertiary spherical aberrations, primary, secondary, tertiary, andquaternary astigmatism. The boundary factor b is

$\begin{matrix}{{b^{2} = {{1740\left( c_{5}^{- 1} \right)^{2}} + \left( c_{5}^{1} \right)^{2} + {1728\left( c_{5}^{- 1} \right)^{2}} + {\left( c_{5}^{1} \right)^{2}\cos \; 2\left( {\theta - \varphi} \right)}}},{where}} & ({A50}) \\{\varphi = {\frac{1}{4}{{\tan^{- 1}\left\lbrack \frac{2c_{5}^{- 1}c_{5}^{1}}{\left( c_{5}^{1} \right)^{2} - \left( c_{5}^{- 1} \right)^{2}} \right\rbrack}.}}} & ({A51})\end{matrix}$

Therefore, propagation of secondary coma becomes elliptical.

H. Quadrafoil Aberration

Finally, for a quadrafoil, Z₄ ⁻⁴ and Z₄ ⁴, the direction factor is

$\begin{matrix}{a = {{{40\left\lbrack {\left( c_{4}^{- 4} \right)^{2} + \left( c_{4}^{4} \right)^{2}} \right\rbrack}Z_{0}^{0}} + {24{\sqrt{3}\left\lbrack {\left( c_{4}^{- 4} \right)^{2} + \left( c_{4}^{4} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {8{\sqrt{5}\left\lbrack {\left( c_{4}^{- 4} \right)^{2} + \left( c_{4}^{4} \right)^{2}} \right\rbrack}Z_{4}^{0}} + {{\frac{8}{\sqrt{7}}\left\lbrack {\left( c_{4}^{- 4} \right)^{2} + \left( c_{4}^{4} \right)^{2}} \right\rbrack}{Z_{6}^{0}.}}}} & ({A52})\end{matrix}$

So the propagation of quadrafoil induces defocus, primary and secondaryspherical aberration. The boundary factor b is

b ²=160[(c ₄ ⁻⁴)²+(c ₄ ²)²].  (A53)

Therefore, propagation of quadrafoil is still circular.

Systems and methods for calculating both the boundary and magnitude of awavefront after it propagates from one plane to another are providedherein. Taylor monomials can be effectively used to achieve theanalytical formulation of the wavefront propagation. Zernikecoefficients can be converted to and from Taylor coefficients forwavefront representation before and after the propagation.

Because of the linear nature of the wavefront as expanded into a set ofbasis functions in some embodiments, the propagation of a wavefront canbe treated using the direction factor and the boundary factor, both ofwhich are not linearly proportional to the wavefront. Therefore, in someembodiments the propagated wavefront is not treated as a linearcombination of the propagation of individual Zernike polynomials. Thepropagation of the low order aberrations can be verified by theclassical vertex correction formula and the propagation of the highorder aberrations can be verified by Zemax® ray tracing software.

The systems and methods disclosed herein are well suited for visioncorrection as the ocular wavefront is measured on one plane and thecorrection is performed on another plane. In some embodiments, theanalytical nature of the results increases the likelihood of a highprecision and in most cases a faster execution.

VIII. Proof of Eq. (A5b) for a Propagated Low Order Wavefront

The low order sphere and cylinder can be expressed in terms of Zernikepolynomials as

$\begin{matrix}\begin{matrix}{{W\left( {{R\; \rho},\theta} \right)} = {{\sqrt{6}c_{2}^{- 2}\rho^{2}\sin \; 2\; \theta} + {\sqrt{3}{c_{2}^{0}\left( {{2\; \rho^{2}} - 1} \right)}} +}} \\{{\sqrt{6}c_{2}^{2}\rho^{2}\cos \; 2\; \theta}} \\{= {{\sqrt{3}{c_{2}^{0}\left( {{2\; \rho^{2}} - 1} \right)}} +}} \\{{{\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}\rho^{2}\cos \; 2\left( {\theta - \varphi} \right)},}}\end{matrix} & \left( {A\; 54} \right)\end{matrix}$

where the cylinder axis φ can be expressed as

$\begin{matrix}{\varphi = {\frac{1}{4}{{\tan^{- 1}\left( \frac{c_{2}^{- 2}}{c_{2}^{2}} \right)}.}}} & ({A55})\end{matrix}$

Without loss of generality, we use a plus cylinder notation in thissection. Therefore, the sphere and cylinder of this wavefront can bederived as

$\begin{matrix}{{S = {{- \frac{4\sqrt{3}c_{2}^{0}}{R^{2}}} - \frac{\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}{R^{2}}}},} & ({A56a}) \\{C = {\frac{4\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}{R^{2}}.}} & ({A56b})\end{matrix}$

Writing Eq. (A54) in Cartesian coordinates, we have

W(u,v)=2√{square root over (6)}c ₂ ⁻² uv+√{square root over (3)}c ₂ ⁰(2u²+2v ²−1)+√{square root over (6)}c ₂ ²(u ² −v ²).  (A57)

Therefore,

$\begin{matrix}\begin{matrix}{{\left\lbrack \frac{\partial{W\left( {u,v} \right)}}{\partial u} \right\rbrack^{2} + \left\lbrack \frac{\partial\left( {u,v} \right)}{\partial v} \right\rbrack^{2}} = {\left( {{2\sqrt{6}c_{2}^{- 2}v} + {4\sqrt{3}c_{2}^{0}u} + {2\sqrt{6}c_{2}^{2}u}} \right)^{2} +}} \\\left. {{2\sqrt{6}c_{2}^{- 2}u} + {4\sqrt{3}c_{2}^{0}v} - {2\sqrt{6}c_{2}^{2}v}} \right)^{2} \\{= {{{24\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}\rho^{2}} +}} \\{{48\sqrt{2}c_{2}^{0}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}\rho^{2}\cos \; 2\left( {\theta - \varphi} \right)}} \\{= {{16\sqrt{3}c_{2}^{0}c_{2}^{- 2}Z_{2}^{- 2}} +}} \\{{{4{\sqrt{3}\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}Z_{2}^{0}} +}} \\{{{16\sqrt{3}c_{2}^{0}c_{2}^{2}Z_{2}^{2}} +}} \\{{{12\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}Z_{0}^{0}}}\end{matrix} & ({A58})\end{matrix}$

Substituting Eq. (A58) into Eq. (A11), we obtain

$\begin{matrix}{{{W^{\prime}\left( {\rho^{\prime},\theta^{\prime}} \right)} = {{{c_{2}^{- 2}Z_{2}^{- 2}} + {c_{2}^{0}Z_{2}^{0}} + {c_{2}^{2}Z_{2}^{2}} - {\frac{d}{2\; R^{2}}\left\{ {{16\sqrt{3}c_{2}^{0}c_{2}^{- 2}Z_{2}^{- 2}} + {4{\sqrt{3}\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}Z_{2}^{0}} + {16\sqrt{3}c_{2}^{0}c_{2}^{2}Z_{2}^{2}} + {{12\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}Z_{0}^{0}}} \right\}}} = {{b_{2}^{- 2}Z_{2}^{- 2}} + {b_{2}^{0}Z_{2}^{0}} + {b_{2}^{2}Z_{2}^{2}} + {b_{0}^{0}Z_{0}^{0}}}}},\mspace{79mu} {where}} & ({A59}) \\{\mspace{79mu} {b_{2}^{- 2} = {\left( {1 - {d\frac{8\sqrt{3}c_{2}^{0}}{R^{2}}}} \right)c_{2}^{- 2}}}} & ({A60a}) \\{\mspace{79mu} {b_{2}^{0} = {\left\{ {1 - {d{\frac{2\sqrt{3}}{c_{2}^{0}R^{2}}\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}}} \right\} c_{2}^{0}}}} & ({A60b}) \\{\mspace{79mu} {b_{2}^{2} = {\left( {1 - {d\frac{8\sqrt{3}c_{2}^{0}}{R^{2}}}} \right)c_{2}^{2}}}} & ({A60c}) \\{\mspace{79mu} {b_{0}^{0} = {{12\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack}.}}} & ({A60d})\end{matrix}$

Hence, the sphere and cylinder of the propagated wavefront are

$\begin{matrix}{{S^{\prime} = {{- \frac{4\sqrt{3}b_{2}^{0}}{R^{\prime 2}}} - \frac{2\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}{R^{\prime 2}}}},} & ({A61a}) \\{C^{\prime} = {\frac{4\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}{R^{\prime 2}}.}} & ({A61b})\end{matrix}$

For the new wavefront radius, we can calculate tan ψ from Eq. (A58)

$\begin{matrix}{{\tan \; \psi} = {\frac{2\sqrt{6}}{R}{\begin{Bmatrix}{\left\lbrack {\left( c_{2}^{- 2} \right)^{2} + {2\left( c_{2}^{0} \right)^{2}} + \left( c_{2}^{2} \right)^{2}} \right\rbrack + {2\sqrt{2}c_{2}^{0}} +} \\{\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}\cos \; 2\left( {\theta - \varphi} \right)}\end{Bmatrix}^{1/2}.}}} & ({A62})\end{matrix}$

Apparently, the shape of the wavefront becomes elliptical from theoriginal circular shape after it propagates a distance d. When θ=φ, theorientation has the minimum power, which corresponds to the spherepower, Eq. (A62) can be written as

$\begin{matrix}{{\tan \; \psi} = {{\frac{1}{R}\left\lbrack {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}.}} & ({A63})\end{matrix}$

Therefore, the semiminor axis of the ellipse is

$\begin{matrix}\begin{matrix}{R_{\min} = {R\left\{ {1 - {\frac{d}{R^{2}}\left\lbrack {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}} \right\}}} \\{= {{R\left( {1 + {dS}} \right)}.}}\end{matrix} & ({A64})\end{matrix}$

Substituting Eq. (A64) and Eq. (A60) into Eq. (A61a) with some algebra,we obtain

$\begin{matrix}\begin{matrix}{S^{\prime} = {{{- \frac{1}{R_{\min}^{2}}}4\sqrt{3}b_{2}^{0}} + {2\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}}} \\{= {- {\frac{1}{R_{\min}^{2}}\begin{bmatrix}{{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}} -} \\{\frac{d}{R^{2}}\left( {{4\sqrt{3}c_{2}^{0}} + {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)}}} \right)^{2}}\end{bmatrix}}}} \\{= {- {\frac{1}{{R^{2}\left( {1 + {dS}} \right)}^{2}}\left\lbrack {{- R^{2}}{S\left( {1 + {dS}} \right)}} \right\rbrack}}} \\{= {\frac{S}{1 + {dS}}.}}\end{matrix} & ({A65})\end{matrix}$

When θ=φ, +π/2, the orientation has the maximum power, which correspondsto combined power of sphere and cylinder, we have

$\begin{matrix}{{\tan \; \psi} = {{\frac{1}{R}\left\lbrack {{4\sqrt{3}c_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}.}} & ({A66})\end{matrix}$

Therefore, the semimajor axis of the ellipse is

$\begin{matrix}\begin{matrix}{R_{\max} = {R\left\{ {1 - {\frac{d}{R^{2}}\left\lbrack {{4\sqrt{3}c_{2}^{0}} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right\rbrack}} \right\}}} \\{= {{R\left\lbrack {1 + {d\left( {S + C} \right)}} \right\rbrack}.}}\end{matrix} & ({A67})\end{matrix}$

Substituting Eq. (A67) and Eq. (A60) into Eq. (A61) with some algebra,we obtain

$\begin{matrix}\begin{matrix}{{S^{\prime} + C^{\prime}} = {{- {\frac{1}{R_{\max}^{2}}\left\lbrack {4\sqrt{3}b_{2}^{0}} \right)}} - {2\sqrt{6}\sqrt{\left( b_{2}^{- 2} \right)^{2} + \left( b_{2}^{2} \right)^{2}}}}} \\{= {- {\frac{1}{R_{\max}^{2}}\begin{bmatrix}{\left. {4\sqrt{3}c_{2}^{0}} \right) - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}} -} \\\left. {{\frac{d}{R^{2}}\left( {4\sqrt{3}c_{2}^{0}} \right)} - {2\sqrt{6}\sqrt{\left( c_{2}^{- 2} \right)^{2} + \left( c_{2}^{2} \right)^{2}}}} \right)^{2}\end{bmatrix}}}} \\{= {- {{\frac{1}{{R^{2}\left\lbrack {1 + {d\left( {S + C} \right)}} \right\rbrack}^{2}}\left\lbrack {{{- R^{2}}S} + C} \right)}\left\lbrack {1 + {d\left( {S + C} \right)}} \right\rbrack}}} \\{= {\frac{S}{1 + {d\left( {S + C} \right)}}.}}\end{matrix} & ({A68})\end{matrix}$

Equations (A65) and (A68) are identical to Eqs. (A5a) and (A5b),respectively, hence proving Eq. (A5b).

IX. Electromagnetic Fields

Embodiments of the present invention also encompass methods and systemsfor evaluating or characterizing an electromagnetic field that ispropagated from a first surface or plane to a second surface or plane.An electromagnetic wave can be represented in a four dimensional space,where one dimensional is the time, and the three other dimensionsrepresent the space. The electric field vector and the magnetic fieldvector are orthogonal to each other, and both of them are orthogonal tothe direction of the propagation of the electromagnetic field. Athree-dimensional field that is defined by the electric field vector andthe magnetic field vector can be written as

$\begin{matrix}{{{\Psi \left( {x,y} \right)} = {{A\left( {x,y} \right)}{\exp \left\lbrack {{- j}\frac{2\; \pi}{\lambda}{\Phi \left( {x,y} \right)}} \right\rbrack}}},} & ({A69})\end{matrix}$

where A(x,y) stands for the modulus, or magnitude, λ is the wavelength,and Φ(x,y) is the phase of the electromagnetic field. In thethree-dimensional space, if the points where the electromagnetic wavehas the same phase are connected, the resulting surface is oftenreferred to as a wavefront. According to some embodiments, if anoriginal plane wave propagates through an isotropic and homogeneousmedium, the electromagnetic wave in a new plane that is parallel to theoriginal plane is in-phase, and consequently there is no wavefronterror. However, if a plane wave propagates through a lens or otheroptical or ocular system, the electromagnetic wave in a new plane thatis parallel to or corresponds to the original plane is typically nolonger in-phase. The difference in terms of the optical path (OPD) inthe three-dimensional space can define a wavefront error.

The energy, or the strength, of the electromagnetic wave may change whenthe wave propagates through a certain medium. Depending upon theproperties of the medium, either or both the magnitude A(x,y) and phaseΦ(x,y) may change. For example, in astronomy, both scintillation andphase fluctuation occur when a plane wave passes through the atmosphericturbulence as both A(x,y) and Φ(x,y) change. For ocular aberrations,A(x,y) in general does not change, or the change can be negligible.Therefore, ocular aberrations are dominated by phase error. For oldereyes, A(x,y) may change significantly, as scattering occurs in thecrystalline lens.

Thus, embodiments of the present invention can involve thecharacterization or evaluation of the phase of an electromagnetic field,or in other words a wavefront. Similarly, the magnitude of theelectromagnetic field is related to scattering in vision application. Insome embodiments, the combination of the magnitude and the phase of thefield can be referred to as an electromagnetic strength, or as astrength of an electromagnetic field. The phase aspect of theelectromagnetic strength involves what is often referred to as awavefront.

As noted above, A(x,y) can represent the magnitude of a complexelectromagnetic field (or wave). In vision analysis, this term is oftenignored or considered to be negligible. That is, A(x,y) is considered asa constant of space and time that does not change as the wavepropagates. On the other hand, Φ(x,y), which is the phase of the complexelectromagnetic wave, typically defines the wavefront, and can also bereferred to as W(x,y). When a wavefront propagates, or in other wordswhen an electromagnetic field or wave propagates, both A(x,y) and Φ(x,y)(or W(x,y)) can be expected to change. Embodiments of the presentinvention encompass methods and systems that can evaluate orcharacterize how Φ(x,y) or W(x,y) changes with a propagation distance(d).

Because Φ(x,y) can represent the wavefront aspect of an electromagneticwave, and can also be referred to as W(x,y), Eq. (A69) can also bewritten as:

$\begin{matrix}{{{\Psi \left( {x,y} \right)} = {{A\left( {x,y} \right)}{\exp \left\lbrack {{- j}\frac{2\; \pi}{\lambda}{W\left( {x,y} \right)}} \right\rbrack}}},} & ({A70})\end{matrix}$

where A(x,y) stands for the modulus, or magnitude, λ is the wavelength,and W(x,y) is the phase of the electromagnetic field. The energy, or thestrength, of the electromagnetic wave may change when the wavepropagates through a certain medium. Depending upon the properties ofthe medium, either or both the magnitude A(x,y) and phase W(x,y) maychange. For example, in astronomy, both scintillation and phasefluctuation occur when a plane wave passes through the atmosphericturbulence as both A(x,y) and W(x,y) change. The wavefront W(x,y) canchange as a function of the propagation distance d. In addition, theboundary of W(x,y) can also changes as it propagates. Similarly, themagnitude of the phase Φ(x,y) or W(x,y) can change as a result ofpropagation. The magnitude of Φ(x,y) or W(x,y) can refer to thenumerical values of Φ(x,y) or W(x,y) at each point of a 2-dimensionalgrid.

A(x,y) can refer to the numerical values of A(x,y) at each point of a2-dimensional grid. As an illustration of one exemplary embodiment, itis helpful to consider a plane wave which can propagate from the retinatowards the pupil of a human eye. If the pupil is sampled with 100×100discrete points, the electromagnetic field at each point is a complexnumber, the modulus A(i,j) with the phase W(i,j), where i and j standfor the indices in this 2-D matrix. If there are no ocular aberrations,the values A and W would be constant. If there are ocular aberrationsbut no scattering, then A can still remain constant but W can bedifferent for different pairs of (i,j). For example, W(i,j) can changeas follows. First, the magnitude (numerical values) can change as thepropagation distance d changes. Second, the boundary can change, or inother words the original 2-dimensional grid can change. The boundarychange can be somewhat like a distortion.

Thus, embodiments of the present invention which encompass methods andsystems for evaluating or characterizing an electromagnetic field thatis propagated from a first surface or plane to a second surface or planecan involve determining a first surface characterization of theelectromagnetic field corresponding to the first surface. The firstsurface characterization can include a first surface field strength.These techniques can also involve determining a propagation distancebetween the first surface and a second surface, and determining a secondsurface characterization of the electromagnetic field based on the firstsurface characterization and the propagation distance, where the secondsurface characterization includes a second surface field strength. Insome cases, the first surface field strength includes a first surfacefield phase, and the second surface field strength includes a secondsurface field phase.

Each of the calculations or operations disclosed herein may be performedusing a computer or other processor having hardware, software, and/orfirmware. The various method steps may be performed by modules, and themodules may comprise any of a wide variety of digital and/or analog dataprocessing hardware and/or software arranged to perform the method stepsdescribed herein. The modules optionally comprising data processinghardware adapted to perform one or more of these steps by havingappropriate machine programming code associated therewith, the modulesfor two or more steps (or portions of two or more steps) beingintegrated into a single processor board or separated into differentprocessor boards in any of a wide variety of integrated and/ordistributed processing architectures. These methods and systems willoften employ a tangible media embodying machine-readable code withinstructions for performing the method steps described above. Suitabletangible media may comprise a memory (including a volatile memory and/ora non-volatile memory), a storage media (such as a magnetic recording ona floppy disk, a hard disk, a tape, or the like; on an optical memorysuch as a CD, a CD-R/W, a CD-ROM, a DVD, or the like; or any otherdigital or analog storage media), or the like.

All patent filings, scientific journals, books, treatises, and otherpublications and materials discussed in this application are herebyincorporated by reference for all purposes. Although embodiments of theinvention have often been described herein with specific reference to awavefront system using lenslets, other suitable wavefront systems thatmeasure angles of light passing through the eye may be employed. Forexample, systems using the principles of ray tracing aberrometry,tscherning aberrometry, and dynamic skiascopy may be used withembodiments of the current invention. The above systems are availablefrom TRACEY Technologies of Bellaire, Tex., Wavelight of Erlangen,Germany, and Nidek, Inc. of Fremont, Calif., respectively. The inventionmay also be practiced with a spatially resolved refractometer asdescribed in U.S. Pat. Nos. 6,099,125; 6,000,800; and 5,258,791, thefull disclosures of which are incorporated herein by reference.Treatments that may benefit from the invention include intraocularlenses, contact lenses, spectacles and other surgical methods inaddition to refractive laser corneal surgery.

While the exemplary embodiments have been described in some detail, byway of example and for clarity of understanding, those of skill in theart will recognize that a variety of modification, adaptations, andchanges may be employed. Hence, the scope of the present inventionshould be limited solely by the appending claims.

1. A method of characterizing an electromagnetic field that ispropagated from a first surface to a second surface, comprising:determining a first surface characterization of the electromagneticfield corresponding to the first surface, the first surfacecharacterization comprising a first surface field strength; determininga propagation distance between the first surface and a second surface;and determining a second surface characterization of the electromagneticfield based on the first surface characterization and the propagationdistance, the second surface characterization comprising a secondsurface field strength.
 2. The method of claim 1, wherein the firstsurface field strength comprises a first surface field phase, and thesecond surface field strength comprises a second surface field phase.